atom_finding.py 38.5 KB
Newer Older
1
2
3
4
# -*- coding: utf-8 -*-
"""
@author: Oleg Ovchinnikov
"""
5
from __future__ import division, print_function, absolute_import
6

7

Unknown's avatar
Unknown committed
8
def apply_select_channel(file_in_h5, img_num, channel_num):
9
10
11
    import numpy as np
    import h5py as h5

Unknown's avatar
Unknown committed
12
13
    main_h5_handle = h5.File(file_in_h5, 'r+')
    image_path = "/Frame_%04i/Channel_%02i" % (img_num, channel_num)
Unknown's avatar
Unknown committed
14
    image_path = "%s/Raw_Data" % image_path
Unknown's avatar
Unknown committed
15
16

    h5_image = main_h5_handle.get(image_path)
17
18
    img2 = np.empty(h5_image.shape, dtype=h5_image.dtype)
    h5_image.read_direct(img2)
Unknown's avatar
Unknown committed
19

Unknown's avatar
Unknown committed
20
21
22
23
24
25
26
27
    sec_i_ref = h5_image.attrs.get("Spectroscopic_Indices")
    sec_i_reg = h5_image.attrs.get("Spectroscopic_Indices_Region")
    sec_v_ref = h5_image.attrs.get("Spectroscopic_Values")
    sec_v_reg = h5_image.attrs.get("Spectroscopic_Values_Region")
    pos_i_ref = h5_image.attrs.get("Position_Indices")
    pos_i_reg = h5_image.attrs.get("Position_Indices_Region")
    pos_v_ref = h5_image.attrs.get("Position_Values")
    pos_v_reg = h5_image.attrs.get("Position_Values_Region")
Unknown's avatar
Unknown committed
28
29
30
    current_ref = h5_image.ref
    current_reg = h5_image.regionref[0:len(img2)]

Unknown's avatar
Unknown committed
31
32
    spec_ind = main_h5_handle[sec_i_ref]
    spec_ind = spec_ind[sec_i_reg]
Unknown's avatar
Unknown committed
33

Unknown's avatar
Unknown committed
34
35
    posi_ind = main_h5_handle[pos_i_ref]
    posi_ind = posi_ind[pos_i_reg]
Unknown's avatar
Unknown committed
36

Unknown's avatar
Unknown committed
37
    image_path = "/Frame_%04i/Channel_Current/Filter_Step_0000" % img_num
Unknown's avatar
Unknown committed
38

39
40
41
    try:
        main_h5_handle.__delitem__(image_path)
    except:
Unknown's avatar
Unknown committed
42
        temp = 1
Unknown's avatar
Unknown committed
43

Unknown's avatar
Unknown committed
44
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
45
46
    main_h5_handle[image_path] = img2
    h5_image_new = main_h5_handle.get(image_path)
47
    h5_new_attrs = h5_image_new.attrs
Unknown's avatar
Unknown committed
48
49
50
51
52
53
54
55
56
57
58
59
60
    h5_new_attrs["Filter_Name"] = "Original"
    h5_new_attrs["Number_Of_Variables"] = 0

    h5_new_attrs["Spectroscopic_Indices"] = sec_i_ref
    h5_new_attrs["Spectroscopic_Indices_Region"] = sec_i_reg
    h5_new_attrs["Spectroscopic_Values"] = sec_v_ref
    h5_new_attrs["Spectroscopic_Values_Region"] = sec_v_reg
    h5_new_attrs["Position_Indices"] = pos_i_ref
    h5_new_attrs["Position_Indices_Region"] = pos_i_reg
    h5_new_attrs["Position_Values"] = pos_v_ref
    h5_new_attrs["Position_Values_Region"] = pos_v_reg
    h5_new_attrs["Parent"] = current_ref
    h5_new_attrs["Parent_Region"] = current_reg
61

Unknown's avatar
Unknown committed
62
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
63
64
65
66

    path_main = main_h5_handle.get(image_path)
    path_attrs = path_main

Unknown's avatar
Unknown committed
67
    channel_name = "Channel_%02i" % channel_num
Unknown's avatar
Unknown committed
68
69
    path_main = main_h5_handle.get(image_path)
    path_attrs = path_main.attrs
Unknown's avatar
Unknown committed
70
71
72
    path_attrs["Origin"] = current_ref
    path_attrs["Origin_Region"] = current_reg
    path_attrs["Origin_Name"] = channel_name
73
74
75
76

    main_h5_handle.close()
    return 1

Unknown's avatar
Unknown committed
77
78

def apply_wiener_filter(file_in_h5, img_num, filter_num):
79
80
    import numpy as np
    import h5py as h5
Unknown's avatar
Unknown committed
81
82

    main_h5_handle = h5.File(file_in_h5, 'r+')
Unknown's avatar
Unknown committed
83
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
84
85
    for x in range(0, filter_num + 1):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)
Unknown's avatar
Unknown committed
86
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
87
88

    h5_image = main_h5_handle.get(image_path)
89
90
    img2 = np.empty(h5_image.shape, dtype=h5_image.dtype)
    h5_image.read_direct(img2)
Unknown's avatar
Unknown committed
91

Unknown's avatar
Unknown committed
92
93
94
95
96
97
98
99
    sec_i_ref = h5_image.attrs.get("Spectroscopic_Indices")
    sec_i_reg = h5_image.attrs.get("Spectroscopic_Indices_Region")
    sec_v_ref = h5_image.attrs.get("Spectroscopic_Values")
    sec_v_reg = h5_image.attrs.get("Spectroscopic_Values_Region")
    pos_i_ref = h5_image.attrs.get("Position_Indices")
    pos_i_reg = h5_image.attrs.get("Position_Indices_Region")
    pos_v_ref = h5_image.attrs.get("Position_Values")
    pos_v_reg = h5_image.attrs.get("Position_Values_Region")
Unknown's avatar
Unknown committed
100
101
102
    current_ref = h5_image.ref
    current_reg = h5_image.regionref[0:len(img2)]

Unknown's avatar
Unknown committed
103
104
    spec_ind = main_h5_handle[sec_i_ref]
    spec_ind = spec_ind[sec_i_reg]
Unknown's avatar
Unknown committed
105

Unknown's avatar
Unknown committed
106
107
    posi_ind = main_h5_handle[pos_i_ref]
    posi_ind = posi_ind[pos_i_reg]
Unknown's avatar
Unknown committed
108
109
110
    img = np.empty([max(posi_ind[:, 0]), max(posi_ind[:, 1])], dtype=h5_image.dtype)
    img[posi_ind[:, 0] - 1, posi_ind[:, 1] - 1] = img2[0:len(img2), 0]

111
    # test shape doe wiener filter
Unknown's avatar
Unknown committed
112
113
114
    s = 1.1  # constant for now change latter

    size = img.shape
115
    h = np.empty(size)
Unknown's avatar
Unknown committed
116
117
118
119
120
121
    width = size[0]
    height = size[1]
    for x in range(0, width):
        for y in range(0, height):
            h[x, y] = (s / 2 * 3.141592) ** (-s * np.sqrt(x ** 2 + y ** 2))
    h = h / np.sum(np.sum(h))
122
    H = np.fft.fft2(h)
Unknown's avatar
Unknown committed
123
    K = np.linspace(.001, 1, 100)
124
    G = np.fft.fft2(img)
Unknown's avatar
Unknown committed
125
126
127
128
129
    errorV = np.empty([1, 100])
    errorV = errorV[0]
    for k1 in range(0, 100):
        W = np.conj(H) / (np.abs(H) ** 2 + K[k1])
        F = W * G
130
        R = np.abs(np.fft.ifft2(F))
Unknown's avatar
Unknown committed
131
132
        errorV[k1] = np.std(np.reshape(R, [1, width * height]))

133
    # use best option 
Unknown's avatar
Unknown committed
134
135
136
137
138
139
140
141
142
    minv = np.min(errorV)
    minl = np.where(errorV == minv)
    minl = minl[0]
    minl = minl[0]
    W = np.conj(H) / (np.abs(H) ** 2 + K[minl])
    F = W * G
    img = np.fft.ifft2(F)

    img = img.reshape(len(img2), 1)
Unknown's avatar
Unknown committed
143
144
    img = abs(img)
    img = np.real(img)
Unknown's avatar
Unknown committed
145

Unknown's avatar
Unknown committed
146
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
147
148
149
    for x in range(0, filter_num + 2):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)

150
151
152
    try:
        main_h5_handle.__delitem__(image_path)
    except:
Unknown's avatar
Unknown committed
153
        temp = 1
Unknown's avatar
Unknown committed
154

Unknown's avatar
Unknown committed
155
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
156
157
    main_h5_handle[image_path] = img
    h5_image_new = main_h5_handle.get(image_path)
158
    h5_new_attrs = h5_image_new.attrs
Unknown's avatar
Unknown committed
159
160
161
162
163
164
165
166
167
168
169
170
171
    h5_new_attrs["Filter_Name"] = "Wiener"
    h5_new_attrs["Number_Of_Variables"] = 0

    h5_new_attrs["Spectroscopic_Indices"] = sec_i_ref
    h5_new_attrs["Spectroscopic_Indices_Region"] = sec_i_reg
    h5_new_attrs["Spectroscopic_Values"] = sec_v_ref
    h5_new_attrs["Spectroscopic_Values_Region"] = sec_v_reg
    h5_new_attrs["Position_Indices"] = pos_i_ref
    h5_new_attrs["Position_Indices_Region"] = pos_i_reg
    h5_new_attrs["Position_Values"] = pos_v_ref
    h5_new_attrs["Position_Values_Region"] = pos_v_reg
    h5_new_attrs["Parent"] = current_ref
    h5_new_attrs["Parent_Region"] = current_reg
Unknown's avatar
Unknown committed
172

173
    main_h5_handle.close()
Unknown's avatar
Unknown committed
174

175
176
    return 1

Unknown's avatar
Unknown committed
177
178

def apply_gaussian_corr_filter(file_in_h5, img_num, filter_num, gauss_width, gauss_box_width):
179
180
    import numpy as np
    import h5py as h5
Unknown's avatar
Unknown committed
181
182

    main_h5_handle = h5.File(file_in_h5, 'r+')
Unknown's avatar
Unknown committed
183
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
184
185
    for x in range(0, filter_num + 1):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)
Unknown's avatar
Unknown committed
186
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
187
188

    h5_image = main_h5_handle.get(image_path)
189
190
    img2 = np.empty(h5_image.shape, dtype=h5_image.dtype)
    h5_image.read_direct(img2)
Unknown's avatar
Unknown committed
191

Unknown's avatar
Unknown committed
192
193
194
195
196
197
198
199
    sec_i_ref = h5_image.attrs.get("Spectroscopic_Indices")
    sec_i_reg = h5_image.attrs.get("Spectroscopic_Indices_Region")
    sec_v_ref = h5_image.attrs.get("Spectroscopic_Values")
    sec_v_reg = h5_image.attrs.get("Spectroscopic_Values_Region")
    pos_i_ref = h5_image.attrs.get("Position_Indices")
    pos_i_reg = h5_image.attrs.get("Position_Indices_Region")
    pos_v_ref = h5_image.attrs.get("Position_Values")
    pos_v_reg = h5_image.attrs.get("Position_Values_Region")
Unknown's avatar
Unknown committed
200
201
202
    current_ref = h5_image.ref
    current_reg = h5_image.regionref[0:len(img2)]

Unknown's avatar
Unknown committed
203
204
    spec_ind = main_h5_handle[sec_i_ref]
    spec_ind = spec_ind[sec_i_reg]
Unknown's avatar
Unknown committed
205

Unknown's avatar
Unknown committed
206
207
    posi_ind = main_h5_handle[pos_i_ref]
    posi_ind = posi_ind[pos_i_reg]
Unknown's avatar
Unknown committed
208
209
210
    img = np.empty([max(posi_ind[:, 0]), max(posi_ind[:, 1])], dtype=h5_image.dtype)
    img[posi_ind[:, 0] - 1, posi_ind[:, 1] - 1] = img2[0:len(img2), 0]

Unknown's avatar
Unknown committed
211
    gauss_cell = [1, 0, 0, gauss_width, gauss_width, 0]
Unknown's avatar
Unknown committed
212

Unknown's avatar
Unknown committed
213
214
215
    s1 = len(img[:, 0])
    s2 = len(img[0, :])
    new_deconv = np.zeros([s1, s2])
Unknown's avatar
Unknown committed
216
217
218

    for k1 in range(0, s1):
        for k2 in range(0, s2):
Unknown's avatar
Unknown committed
219
220
221
222
223
            x_min = max([-gauss_box_width, -k1])
            x_max = min([gauss_box_width, (s1 - k1 - 1)])
            y_min = max([-gauss_box_width, -k2])
            y_max = min([gauss_box_width, (s2 - k2 - 1)])
            [xx, yy] = np.meshgrid(range(y_min, y_max + 1), range(x_min, x_max + 1))
Unknown's avatar
Unknown committed
224
225
226
            gaus = fun_2d_gaussian(xx, yy, gauss_cell)
            temp = np.corrcoef(gaus.reshape([1, gaus.size]),
                               img[k1 + x_min:k1 + x_max + 1,
Unknown's avatar
Unknown committed
227
                               k2 + y_min:k2 + y_max + 1].reshape([1, gaus.size]))
Unknown's avatar
Unknown committed
228
229
            new_deconv[k1, k2] = temp[0, 1]

Unknown's avatar
Unknown committed
230
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
231
232
233
    for x in range(0, filter_num + 2):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)

234
235
236
    try:
        main_h5_handle.__delitem__(image_path)
    except:
Unknown's avatar
Unknown committed
237
        temp = 1
Unknown's avatar
Unknown committed
238

Unknown's avatar
Unknown committed
239
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
240
241
    main_h5_handle[image_path] = new_deconv
    h5_image_new = main_h5_handle.get(image_path)
242
    h5_new_attrs = h5_image_new.attrs
Unknown's avatar
Unknown committed
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
    h5_new_attrs["Filter_Name"] = "2D Correlation"
    h5_new_attrs["Number_Of_Variables"] = 2
    h5_new_attrs["Variable_1_Name"] = "Filter width"
    h5_new_attrs["Variable_1_Value"] = gauss_box_width
    h5_new_attrs["Variable_2_Name"] = "Gassian width"
    h5_new_attrs["Variable_2_Value"] = gauss_width

    h5_new_attrs["Spectroscopic_Indices"] = sec_i_ref
    h5_new_attrs["Spectroscopic_Indices_Region"] = sec_i_reg
    h5_new_attrs["Spectroscopic_Values"] = sec_v_ref
    h5_new_attrs["Spectroscopic_Values_Region"] = sec_v_reg
    h5_new_attrs["Position_Indices"] = pos_i_ref
    h5_new_attrs["Position_Indices_Region"] = pos_i_reg
    h5_new_attrs["Position_Values"] = pos_v_ref
    h5_new_attrs["Position_Values_Region"] = pos_v_reg
    h5_new_attrs["Parent"] = current_ref
    h5_new_attrs["Parent_Region"] = current_reg
Unknown's avatar
Unknown committed
260

261
262
263
    main_h5_handle.close()

    return 1
Unknown's avatar
Unknown committed
264
265
266
267
268


def fun_2d_gaussian(x, y, parm):
    import numpy as np

Unknown's avatar
Unknown committed
269
    amp = np.double(parm[0])
Unknown's avatar
Unknown committed
270

Unknown's avatar
Unknown committed
271
272
    x_cent = np.double(parm[1])
    y_cent = np.double(parm[2])
Unknown's avatar
Unknown committed
273

Unknown's avatar
Unknown committed
274
275
    x_wid = np.double(parm[3])
    y_wid = np.double(parm[4])
Unknown's avatar
Unknown committed
276

Unknown's avatar
Unknown committed
277
    ang = np.double(parm[5])
Unknown's avatar
Unknown committed
278

Unknown's avatar
Unknown committed
279
280
281
    a = ((np.cos(ang) ** 2) / (2 * x_wid ** 2)) + ((np.sin(ang) ** 2) / (2 * y_wid ** 2))
    b = -((np.sin(2 * ang)) / (4 * x_wid ** 2)) + ((np.sin(2 * ang)) / (4 * y_wid ** 2))
    c = ((np.sin(ang) ** 2) / (2 * x_wid ** 2)) + ((np.cos(ang) ** 2) / (2 * y_wid ** 2))
Unknown's avatar
Unknown committed
282
283

    gaussian = amp * (
Unknown's avatar
Unknown committed
284
        np.exp(-((a * (x - x_cent) ** 2) + (2 * b * (x - x_cent) * (y - y_cent)) + (c * (y - y_cent) ** 2))))
285
286
287

    return gaussian

Unknown's avatar
Unknown committed
288
289

def apply_invert_filter(file_in_h5, img_num, filter_num):
290
291
    import numpy as np
    import h5py as h5
Unknown's avatar
Unknown committed
292
293

    main_h5_handle = h5.File(file_in_h5, 'r+')
Unknown's avatar
Unknown committed
294
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
295
296
    for x in range(0, filter_num + 1):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)
Unknown's avatar
Unknown committed
297
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
298
299

    h5_image = main_h5_handle.get(image_path)
300
301
    img = np.empty(h5_image.shape, dtype=h5_image.dtype)
    h5_image.read_direct(img)
Unknown's avatar
Unknown committed
302

Unknown's avatar
Unknown committed
303
304
305
306
307
308
309
310
    sec_i_ref = h5_image.attrs.get("Spectroscopic_Indices")
    sec_i_reg = h5_image.attrs.get("Spectroscopic_Indices_Region")
    sec_v_ref = h5_image.attrs.get("Spectroscopic_Values")
    sec_v_reg = h5_image.attrs.get("Spectroscopic_Values_Region")
    pos_i_ref = h5_image.attrs.get("Position_Indices")
    pos_i_reg = h5_image.attrs.get("Position_Indices_Region")
    pos_v_ref = h5_image.attrs.get("Position_Values")
    pos_v_reg = h5_image.attrs.get("Position_Values_Region")
Unknown's avatar
Unknown committed
311
312
313
    current_ref = h5_image.ref
    current_reg = h5_image.regionref[0:len(img)]

Unknown's avatar
Unknown committed
314
315
316
317
    m_img = img.mean()
    img = img - m_img
    img = -img
    img = img + m_img
Unknown's avatar
Unknown committed
318

Unknown's avatar
Unknown committed
319
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
320
321
322
323
    for x in range(0, filter_num + 2):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)

    temp = 0
324
325
326
    try:
        main_h5_handle.__delitem__(image_path)
    except:
Unknown's avatar
Unknown committed
327
        temp = 1
Unknown's avatar
Unknown committed
328

Unknown's avatar
Unknown committed
329
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
330
331
    main_h5_handle[image_path] = img
    h5_image_new = main_h5_handle.get(image_path)
332
    h5_new_attrs = h5_image_new.attrs
Unknown's avatar
Unknown committed
333
334
335
336
337
338
339
340
341
342
343
344
345
    h5_new_attrs["Filter_Name"] = "invert z contrast"
    h5_new_attrs["Number_Of_Variables"] = 0

    h5_new_attrs["Spectroscopic_Indices"] = sec_i_ref
    h5_new_attrs["Spectroscopic_Indices_Region"] = sec_i_reg
    h5_new_attrs["Spectroscopic_Values"] = sec_v_ref
    h5_new_attrs["Spectroscopic_Values_Region"] = sec_v_reg
    h5_new_attrs["Position_Indices"] = pos_i_ref
    h5_new_attrs["Position_Indices_Region"] = pos_i_reg
    h5_new_attrs["Position_Values"] = pos_v_ref
    h5_new_attrs["Position_Values_Region"] = pos_v_reg
    h5_new_attrs["Parent"] = current_ref
    h5_new_attrs["Parent_Region"] = current_reg
Unknown's avatar
Unknown committed
346

Unknown's avatar
Unknown committed
347
    main_h5_handle.close()
Unknown's avatar
Unknown committed
348

349
    return 1
Unknown's avatar
Unknown committed
350
351
352


def apply_find(file_path_h5, file_name_h5, file_path_png, file_name_png, filter_width, img_num, filter_num):
353
354
    import numpy as np
    import h5py as h5
Unknown's avatar
Unknown committed
355

Unknown's avatar
Unknown committed
356
    image_path = "/Frame_%04i/Filtered_Data/Stack_0000" % img_num
Unknown's avatar
Unknown committed
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
    for x in range(0, filter_num + 1):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)

    s = 1.1  # constant for now change latter

    size = img.shape
    s1 = size[0]
    s2 = size[1]

    mat_large = np.empty([s1 + 2 * filter_width + 1, s2 + 2 * filter_width + 1])
    mat_large[:, :] = np.inf

    for k1 in range(-filter_width, filter_width + 1):
        for k2 in range(-filter_width, filter_width + 1):
            mat_large[filter_width - k1:-(filter_width + k1) - 1,
Unknown's avatar
Unknown committed
372
373
374
375
376
                      filter_width - k2:-(filter_width + k2) - 1] = np.minimum(mat_large[filter_width - k1:
                                                                                         -filter_width - k1 - 1,
                                                                                         filter_width - k2:
                                                                                         -filter_width - k2 - 1],
                                                                                h5_image)
Unknown's avatar
Unknown committed
377
378
379

    deconv_mat_temp = mat_large[filter_width:len(mat_larg[1, :]) - filter_width,
                      filter_width:len(mat_larg[:, 1]) - filter_width]
Unknown's avatar
Unknown committed
380
    filtered_image = h5_image - deconv_mat_temp
381
382

    return 1
Unknown's avatar
Unknown committed
383
384
385


def apply_binarization_filter(file_in_h5, img_num, filter_num):
386
387
    import numpy as np
    import h5py as h5
Unknown's avatar
Unknown committed
388
389

    main_h5_handle = h5.File(file_in_h5, 'r+')
Unknown's avatar
Unknown committed
390
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
391
392
    for x in range(0, filter_num + 1):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)
Unknown's avatar
Unknown committed
393
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
394
395

    h5_image = main_h5_handle.get(image_path)
396
397
    img = np.empty(h5_image.shape, dtype=h5_image.dtype)
    h5_image.read_direct(img)
Unknown's avatar
Unknown committed
398

Unknown's avatar
Unknown committed
399
400
401
402
403
404
405
406
    sec_i_ref = h5_image.attrs.get("Spectroscopic_Indices")
    sec_i_reg = h5_image.attrs.get("Spectroscopic_Indices_Region")
    sec_v_ref = h5_image.attrs.get("Spectroscopic_Values")
    sec_v_reg = h5_image.attrs.get("Spectroscopic_Values_Region")
    pos_i_ref = h5_image.attrs.get("Position_Indices")
    pos_i_reg = h5_image.attrs.get("Position_Indices_Region")
    pos_v_ref = h5_image.attrs.get("Position_Values")
    pos_v_reg = h5_image.attrs.get("Position_Values_Region")
Unknown's avatar
Unknown committed
407
408
409
    current_ref = h5_image.ref
    current_reg = h5_image.regionref[0:len(img)]

Unknown's avatar
Unknown committed
410
411
    spec_ind = main_h5_handle[sec_i_ref]
    spec_ind = spec_ind[sec_i_reg]
Unknown's avatar
Unknown committed
412

Unknown's avatar
Unknown committed
413
414
    posi_ind = main_h5_handle[pos_i_ref]
    posi_ind = posi_ind[pos_i_reg]
Unknown's avatar
Unknown committed
415
416
417
418
419
420
421
422
423
424
425
426
427
428

    i_max = img.max()
    i_min = img.min()
    i_diff = i_max - i_min

    max_r = 49
    filter_img = np.zeros([max_r + 1, len(img)])
    time_out = np.zeros([max_r + 1, 1])
    time_out_i = np.zeros([max_r + 1, 1])

    for x in range(0, max_r + 1):
        temp = np.zeros([len(img), 1])
        r = x / max_r
        temp[img > (i_min + (i_diff * r))] = 1
Unknown's avatar
Unknown committed
429
430
431
        filter_img[x, :] = temp[:, 0]
        time_out[x] = r
        time_out_i[x] = x + 1
Unknown's avatar
Unknown committed
432

Unknown's avatar
Unknown committed
433
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
434
435
436
    for x in range(0, filter_num + 2):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)

437
438
439
    try:
        main_h5_handle.__delitem__(image_path)
    except:
Unknown's avatar
Unknown committed
440
        temp = 1
Unknown's avatar
Unknown committed
441

Unknown's avatar
Unknown committed
442
    image_path_f = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
443
444
    main_h5_handle[image_path_f] = img
    h5_image_new = main_h5_handle.get(image_path_f)
445
    h5_new_attrs = h5_image_new.attrs
Unknown's avatar
Unknown committed
446
447
448
449
450
451
452
453
454
455
456
457
458
    h5_new_attrs["Filter_Name"] = "Threshold_Reference_Image"
    h5_new_attrs["Number_Of_Variables"] = 0

    h5_new_attrs["Spectroscopic_Indices"] = sec_i_ref
    h5_new_attrs["Spectroscopic_Indices_Region"] = sec_i_reg
    h5_new_attrs["Spectroscopic_Values"] = sec_v_ref
    h5_new_attrs["Spectroscopic_Values_Region"] = sec_v_reg
    h5_new_attrs["Position_Indices"] = pos_i_ref
    h5_new_attrs["Position_Indices_Region"] = pos_i_reg
    h5_new_attrs["Position_Values"] = pos_v_ref
    h5_new_attrs["Position_Values_Region"] = pos_v_reg
    h5_new_attrs["Parent"] = current_ref
    h5_new_attrs["Parent_Region"] = current_reg
Unknown's avatar
Unknown committed
459

Unknown's avatar
Unknown committed
460
    image_path_sv = "%s/Spectroscopic_Values" % image_path
Unknown's avatar
Unknown committed
461
462
    main_h5_handle[image_path_sv] = time_out
    h5_image_new = main_h5_handle.get(image_path_sv)
Unknown's avatar
Unknown committed
463
464
465
    new_sv_ref = h5_image_new.ref
    new_sv_reg = h5_image_new.regionref[0:len(time_out)]

Unknown's avatar
Unknown committed
466
    image_path_si = "%s/Spectroscopic_Indices" % image_path
Unknown's avatar
Unknown committed
467
468
    main_h5_handle[image_path_si] = time_out_i
    h5_image_new = main_h5_handle.get(image_path_si)
Unknown's avatar
Unknown committed
469
470
471
    new_si_ref = h5_image_new.ref
    new_si_reg = h5_image_new.regionref[0:len(time_out)]

Unknown's avatar
Unknown committed
472
    image_path_b = "%s/Binary_Matrix" % image_path
Unknown's avatar
Unknown committed
473
474
    main_h5_handle[image_path_b] = filter_img
    h5_image_new = main_h5_handle.get(image_path_b)
475
    h5_new_attrs = h5_image_new.attrs
Unknown's avatar
Unknown committed
476
477
478
479
480
481
482
483
484
485
486
487
488
    h5_new_attrs["Filter_Name"] = "Threshold_Matrix"
    h5_new_attrs["Number_Of_Variables"] = 0

    h5_new_attrs["Spectroscopic_Indices"] = new_si_ref
    h5_new_attrs["Spectroscopic_Indices_Region"] = new_si_reg
    h5_new_attrs["Spectroscopic_Values"] = new_sv_ref
    h5_new_attrs["Spectroscopic_Values_Region"] = new_sv_reg
    h5_new_attrs["Position_Indices"] = pos_i_ref
    h5_new_attrs["Position_Indices_Region"] = pos_i_reg
    h5_new_attrs["Position_Values"] = pos_v_ref
    h5_new_attrs["Position_Values_Region"] = pos_v_reg
    h5_new_attrs["Parent"] = current_ref
    h5_new_attrs["Parent_Region"] = current_reg
489
490
491
492

    main_h5_handle.close()

    return 1
Unknown's avatar
Unknown committed
493
494
495


def apply_binarization_filter_select(file_in_h5, img_num, filter_num, threshold):
496
497
    import numpy as np
    import h5py as h5
Unknown's avatar
Unknown committed
498
499

    main_h5_handle = h5.File(file_in_h5, 'r+')
Unknown's avatar
Unknown committed
500
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
501
502
    for x in range(0, filter_num + 1):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)
Unknown's avatar
Unknown committed
503
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
504
505

    h5_image = main_h5_handle.get(image_path)
506
507
    img = np.empty(h5_image.shape, dtype=h5_image.dtype)
    h5_image.read_direct(img)
Unknown's avatar
Unknown committed
508

Unknown's avatar
Unknown committed
509
510
511
512
513
514
515
516
    sec_i_ref = h5_image.attrs.get("Spectroscopic_Indices")
    sec_i_reg = h5_image.attrs.get("Spectroscopic_Indices_Region")
    sec_v_ref = h5_image.attrs.get("Spectroscopic_Values")
    sec_v_reg = h5_image.attrs.get("Spectroscopic_Values_Region")
    pos_i_ref = h5_image.attrs.get("Position_Indices")
    pos_i_reg = h5_image.attrs.get("Position_Indices_Region")
    pos_v_ref = h5_image.attrs.get("Position_Values")
    pos_v_reg = h5_image.attrs.get("Position_Values_Region")
Unknown's avatar
Unknown committed
517
518
519
    current_ref = h5_image.ref
    current_reg = h5_image.regionref[0:len(img)]

Unknown's avatar
Unknown committed
520
521
    spec_ind = main_h5_handle[sec_i_ref]
    spec_ind = spec_ind[sec_i_reg]
Unknown's avatar
Unknown committed
522

Unknown's avatar
Unknown committed
523
524
    posi_ind = main_h5_handle[pos_i_ref]
    posi_ind = posi_ind[pos_i_reg]
Unknown's avatar
Unknown committed
525
526
527
528
529
530
531
532
533
534

    i_max = img.max()
    i_min = img.min()
    i_diff = i_max - i_min

    max_r = 49.0
    filter_img = np.zeros([1, len(img)])
    temp = np.zeros([len(img), 1])
    r = threshold / max_r
    temp[img > (i_min + (i_diff * r))] = 1
Unknown's avatar
Unknown committed
535
    filter_img = temp[:, 0]
Unknown's avatar
Unknown committed
536

Unknown's avatar
Unknown committed
537
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
538
539
540
    for x in range(0, filter_num + 2):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)

541
542
543
    try:
        main_h5_handle.__delitem__(image_path)
    except:
Unknown's avatar
Unknown committed
544
        temp = 1
Unknown's avatar
Unknown committed
545

Unknown's avatar
Unknown committed
546
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
547
548
    main_h5_handle[image_path] = filter_img
    h5_image_new = main_h5_handle.get(image_path)
549
    h5_new_attrs = h5_image_new.attrs
Unknown's avatar
Unknown committed
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
    h5_new_attrs["Filter_Name"] = "Threshold"
    h5_new_attrs["Number_Of_Variables"] = 1
    h5_new_attrs["Variable_1_Name"] = "Threshold"
    h5_new_attrs["Variable_1_Value"] = r

    h5_new_attrs["Spectroscopic_Indices"] = sec_i_ref
    h5_new_attrs["Spectroscopic_Indices_Region"] = sec_i_reg
    h5_new_attrs["Spectroscopic_Values"] = sec_v_ref
    h5_new_attrs["Spectroscopic_Values_Region"] = sec_v_reg
    h5_new_attrs["Position_Indices"] = pos_i_ref
    h5_new_attrs["Position_Indices_Region"] = pos_i_reg
    h5_new_attrs["Position_Values"] = pos_v_ref
    h5_new_attrs["Position_Values_Region"] = pos_v_reg
    h5_new_attrs["Parent"] = current_ref
    h5_new_attrs["Parent_Region"] = current_reg
565
566
567
568
569
570

    main_h5_handle.close()

    return


Unknown's avatar
Unknown committed
571
572
def cluster_into_atomic_columns(file_in_h5, img_num, filter_num, dist_val):
    import numpy as np
573
    import h5py as h5
Unknown's avatar
Unknown committed
574
575

    main_h5_handle = h5.File(file_in_h5, 'r+')
Unknown's avatar
Unknown committed
576
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
577
578
    for ifilt in range(0, filter_num + 1):
        image_path = "%s/Filter_Step_%04i" % (image_path, ifilt)
Unknown's avatar
Unknown committed
579
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
580
581

    h5_image = main_h5_handle.get(image_path)
582
583
    img2 = np.empty(h5_image.shape, dtype=h5_image.dtype)
    h5_image.read_direct(img2)
Unknown's avatar
Unknown committed
584

Unknown's avatar
Unknown committed
585
586
587
588
589
590
591
592
    sec_i_ref = h5_image.attrs.get("Spectroscopic_Indices")
    sec_i_reg = h5_image.attrs.get("Spectroscopic_Indices_Region")
    sec_v_ref = h5_image.attrs.get("Spectroscopic_Values")
    sec_v_reg = h5_image.attrs.get("Spectroscopic_Values_Region")
    pos_i_ref = h5_image.attrs.get("Position_Indices")
    pos_i_reg = h5_image.attrs.get("Position_Indices_Region")
    pos_v_ref = h5_image.attrs.get("Position_Values")
    pos_v_reg = h5_image.attrs.get("Position_Values_Region")
Unknown's avatar
Unknown committed
593
594
595
    current_ref_m = h5_image.ref
    current_reg_m = h5_image.regionref[0:len(img2)]

Unknown's avatar
Unknown committed
596
597
    spec_ind = main_h5_handle[sec_i_ref]
    spec_ind = spec_ind[sec_i_reg]
Unknown's avatar
Unknown committed
598

Unknown's avatar
Unknown committed
599
600
    posi_ind = main_h5_handle[pos_i_ref]
    posi_ind = posi_ind[pos_i_reg]
Unknown's avatar
Unknown committed
601
    img = np.empty([max(posi_ind[:, 0]), max(posi_ind[:, 1])], dtype=h5_image.dtype)
602
    try:
Unknown's avatar
Unknown committed
603
        img[posi_ind[:, 0] - 1, posi_ind[:, 1] - 1] = img2[0:len(img2), 0]
604
    except:
Unknown's avatar
Unknown committed
605
        img[posi_ind[:, 0] - 1, posi_ind[:, 1] - 1] = img2[0:len(img2)]
606

Unknown's avatar
Unknown committed
607
608
    centers = cluster_2d_oleg_return_geo_center(img, dist_val)

Unknown's avatar
Unknown committed
609
610
    image_path_org = "/Frame_%04i/Channel_Current" % img_num
    image_path_new = "/Frame_%04i/Channel_Finished" % img_num
611
612
613
614

    try:
        main_h5_handle.__delitem__(image_path_new)
    except:
Unknown's avatar
Unknown committed
615
        temp = 1
Unknown's avatar
Unknown committed
616
617

    old_folder = main_h5_handle.get(image_path_org)
618
    main_h5_handle.create_group(image_path_new)
Unknown's avatar
Unknown committed
619
    new_folder = main_h5_handle.get(image_path_new)
620
621
    new_folder_attrs = new_folder.attrs

Unknown's avatar
Unknown committed
622
623
624
625
626
627
    parrent_ref = old_folder.attrs.get("Origin")
    parrent_reg = old_folder.attrs.get("Origin_Region")
    parrent_name = old_folder.attrs.get("Origin_Name")
    new_folder_attrs["Origin"] = parrent_ref
    new_folder_attrs["Origin_Region"] = parrent_reg
    new_folder_attrs["Origin_Name"] = parrent_name
Unknown's avatar
Unknown committed
628
629
630
631

    for ifilt in range(0, filter_num + 1):
        image_path_org = "%s/Filter_Step_%04i" % (image_path_org, ifilt)
        image_path_new = "%s/Filter_Step_%04i" % (image_path_new, ifilt)
Unknown's avatar
Unknown committed
632
633
        image_path_org_temp = "%s/Filtered_Image" % image_path_org
        image_path_new_temp = "%s/Filtered_Image" % image_path_new
Unknown's avatar
Unknown committed
634
635

        h5_image_old = main_h5_handle.get(image_path_org_temp)
636
637
        img_old = np.empty(h5_image_old.shape, dtype=h5_image_old.dtype)
        h5_image_old.read_direct(img_old)
Unknown's avatar
Unknown committed
638

Unknown's avatar
Unknown committed
639
640
641
642
643
644
645
646
647
648
649
650
651
        sec_i_ref = h5_image_old.attrs.get("Spectroscopic_Indices")
        sec_i_reg = h5_image_old.attrs.get("Spectroscopic_Indices_Region")
        sec_v_ref = h5_image_old.attrs.get("Spectroscopic_Values")
        sec_v_reg = h5_image_old.attrs.get("Spectroscopic_Values_Region")
        pos_i_ref = h5_image_old.attrs.get("Position_Indices")
        pos_i_reg = h5_image_old.attrs.get("Position_Indices_Region")
        pos_v_ref = h5_image_old.attrs.get("Position_Values")
        pos_v_rexg = h5_image_old.attrs.get("Position_Values_Region")
        filter_name = h5_image_old.attrs.get("Filter_Name")
        number_var = h5_image_old.attrs.get("Number_Of_Variables")

        main_h5_handle[image_path_new_temp] = img_old
        h5_image_new = main_h5_handle.get(image_path)
652
        h5_new_attrs = h5_image_new.attrs
Unknown's avatar
Unknown committed
653
654
        h5_new_attrs["Filter_Name"] = filter_name
        h5_new_attrs["Number_Of_Variables"] = number_var
Unknown's avatar
Unknown committed
655
656

        for ivar in range(1, number_var + 1):
Unknown's avatar
Unknown committed
657
658
659
660
            var_name = h5_image_old.attrs.get("Variable_%01i_Name" % ivar)
            var_value = h5_image_old.attrs.get("Variable_%01i_Value" % ivar)
            h5_new_attrs["Variable_%01i_Name" % ivar] = var_name
            h5_new_attrs["Variable_%01i_Value" % ivar] = var_value
Unknown's avatar
Unknown committed
661
662
663
664
665
666
667
668
669
670
671

        h5_new_attrs["Spectroscopic_Indices"] = sec_i_ref
        h5_new_attrs["Spectroscopic_Indices_Region"] = sec_i_reg
        h5_new_attrs["Spectroscopic_Values"] = sec_v_ref
        h5_new_attrs["Spectroscopic_Values_Region"] = sec_v_reg
        h5_new_attrs["Position_Indices"] = pos_i_ref
        h5_new_attrs["Position_Indices_Region"] = pos_i_reg
        h5_new_attrs["Position_Values"] = pos_v_ref
        h5_new_attrs["Position_Values_Region"] = pos_v_reg
        h5_new_attrs["Parent"] = parrent_ref
        h5_new_attrs["Parent_Region"] = parrent_reg
Unknown's avatar
Unknown committed
672
673

        parrent_ref = h5_image_new.ref
Unknown's avatar
Unknown committed
674
        parrent_reg = h5_image_new.regionref[0:len(img_old)]
Unknown's avatar
Unknown committed
675

Unknown's avatar
Unknown committed
676
    image_path = "%s/Lattice/Positions" % image_path_new
Unknown's avatar
Unknown committed
677
678
    main_h5_handle[image_path] = centers
    h5_image_new = main_h5_handle.get(image_path)
679
    h5_new_attrs = h5_image_new.attrs
Unknown's avatar
Unknown committed
680
681
682
683
684
685
    h5_new_attrs["Method"] = "Density_Clustering"
    h5_new_attrs["Number_Of_Variables"] = 1
    h5_new_attrs["Variable_1_Name"] = "Connection_Distance"
    h5_new_attrs["Variable_1_Value"] = dist_val
    h5_new_attrs["Parent"] = current_ref_m
    h5_new_attrs["Parent_Region"] = current_reg_m
Unknown's avatar
Unknown committed
686

687
688
    main_h5_handle.close()

Unknown's avatar
Unknown committed
689
690

def cluster_2d_oleg(mat_in, dist_val):
691
    import numpy as np
Unknown's avatar
Unknown committed
692
693

    tt = [0, 0]
Unknown's avatar
Unknown committed
694
695
    clusters = []
    to_cluster = np.argwhere(mat_in)
Unknown's avatar
Unknown committed
696
697
    to_cluster = to_cluster.tolist()

Unknown's avatar
Unknown committed
698
    while len(to_cluster) > 0:
Unknown's avatar
Unknown committed
699
700
701
        clust = []
        final_clust = []
        clust.append(to_cluster[0])
702
        to_cluster.remove(to_cluster[0])
Unknown's avatar
Unknown committed
703
        while (len(clust) > 0) & (len(to_cluster) > 0):
Unknown's avatar
Unknown committed
704
705
706
707
708
709
710
            tt[0] = 5000.0 * dist_val
            tt[1] = len(to_cluster)
            t1 = min(tt)
            t1 = int(t1)
            to_cluster_t = to_cluster[0:t1]
            dem_diff = abs(to_cluster_t - np.tile(clust[0], [t1, 1]))
            diff_vec = np.argwhere(sum(dem_diff, 1) <= dist_val)
Unknown's avatar
Unknown committed
711
            for px in range(len(diff_vec), 0, -1):
Unknown's avatar
Unknown committed
712
713
                pt = diff_vec[px - 1]
                clust.append(to_cluster[pt])
714
                to_cluster.remove(to_cluster[pt])
Unknown's avatar
Unknown committed
715
716
            final_clust.append(clust[0])
            clust.remove(clust[0])
Unknown's avatar
Unknown committed
717

Unknown's avatar
Unknown committed
718
719
        if len(clust) > 0:
            while len(clust) < 0:
Unknown's avatar
Unknown committed
720
721
                final_clust.append(clust[0])
                clust.remove(clust[0])
Unknown's avatar
Unknown committed
722

Unknown's avatar
Unknown committed
723
        clusters.append(final_clust)
Unknown's avatar
Unknown committed
724

725
726
    return clusters

Unknown's avatar
Unknown committed
727
728

def cluster_2d_oleg_return_geo_center(mat_in, dist_val):
729
    import numpy as np
Unknown's avatar
Unknown committed
730
731

    tt = [0, 0]
Unknown's avatar
Unknown committed
732
733
    clusters = []
    to_cluster = np.argwhere(mat_in)
Unknown's avatar
Unknown committed
734
735
    to_cluster = to_cluster.tolist()

Unknown's avatar
Unknown committed
736
    while len(to_cluster) > 0:
Unknown's avatar
Unknown committed
737
738
739
        clust = []
        final_clust = []
        clust.append(to_cluster[0])
740
        to_cluster.remove(to_cluster[0])
Unknown's avatar
Unknown committed
741
        while (len(clust) > 0) & (len(to_cluster) > 0):
Unknown's avatar
Unknown committed
742
743
744
745
746
747
748
            tt[0] = 5000.0 * dist_val
            tt[1] = len(to_cluster)
            t1 = min(tt)
            t1 = int(t1)
            to_cluster_t = to_cluster[0:t1]
            dem_diff = abs(to_cluster_t - np.tile(clust[0], [t1, 1]))
            diff_vec = np.argwhere(np.sum(dem_diff, axis=1) <= dist_val)
Unknown's avatar
Unknown committed
749
            for px in range(len(diff_vec), 0, -1):
Unknown's avatar
Unknown committed
750
751
                pt = diff_vec[px - 1]
                clust.append(to_cluster[pt])
752
                to_cluster.remove(to_cluster[pt])
Unknown's avatar
Unknown committed
753
754
            final_clust.append(clust[0])
            clust.remove(clust[0])
Unknown's avatar
Unknown committed
755

Unknown's avatar
Unknown committed
756
757
        if len(clust) > 0:
            while len(clust) < 0:
Unknown's avatar
Unknown committed
758
759
                final_clust.append(clust[0])
                clust.remove(clust[0])
760

Unknown's avatar
Unknown committed
761
        if len(final_clust) > 2:
Unknown's avatar
Unknown committed
762
            clusters.append(np.mean(final_clust, 0))
Unknown's avatar
Unknown committed
763
764

    clusters = np.array(clusters)
Unknown's avatar
Unknown committed
765
    clusters = clusters.tolist()
Unknown's avatar
Unknown committed
766

767
768
769
    return clusters


Unknown's avatar
Unknown committed
770
def return_img(file_in_h5, img_num, filter_num):
771
772
    import numpy as np
    import h5py as h5
Unknown's avatar
Unknown committed
773
774

    main_h5_handle = h5.File(file_in_h5, 'r+')
Unknown's avatar
Unknown committed
775
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
776
777
    for x in range(0, filter_num + 1):
        image_path = "%s/Filter_Step_%04i" % (image_path, x)
Unknown's avatar
Unknown committed
778
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
779
780

    h5_image = main_h5_handle.get(image_path)
781
782
    img2 = np.empty(h5_image.shape, dtype=h5_image.dtype)
    h5_image.read_direct(img2)
Unknown's avatar
Unknown committed
783

Unknown's avatar
Unknown committed
784
785
786
787
788
789
790
791
    sec_i_ref = h5_image.attrs.get("Spectroscopic_Indices")
    sec_i_reg = h5_image.attrs.get("Spectroscopic_Indices_Region")
    sec_v_ref = h5_image.attrs.get("Spectroscopic_Values")
    sec_v_reg = h5_image.attrs.get("Spectroscopic_Values_Region")
    pos_i_ref = h5_image.attrs.get("Position_Indices")
    pos_i_reg = h5_image.attrs.get("Position_Indices_Region")
    pos_v_ref = h5_image.attrs.get("Position_Values")
    pos_v_reg = h5_image.attrs.get("Position_Values_Region")
Unknown's avatar
Unknown committed
792

Unknown's avatar
Unknown committed
793
794
    spec_ind = main_h5_handle[sec_i_ref]
    spec_ind = spec_ind[sec_i_reg]
Unknown's avatar
Unknown committed
795

Unknown's avatar
Unknown committed
796
797
    posi_ind = main_h5_handle[pos_i_ref]
    posi_ind = posi_ind[pos_i_reg]
798
799
    img = np.empty([max(posi_ind[:, 0]), max(posi_ind[:, 1])], dtype=h5_image.dtype)
    try:
Unknown's avatar
Unknown committed
800
        img[posi_ind[:, 0] - 1, posi_ind[:, 1] - 1] = img2[0:len(img2), 0]
801
    except:
Unknown's avatar
Unknown committed
802
        img[posi_ind[:, 0] - 1, posi_ind[:, 1] - 1] = img2[0:len(img2)]
803
804
805
    main_h5_handle.close()

    return img
Unknown's avatar
Unknown committed
806
807
808


def return_pos(file_in_h5, img_num):
809
810
    import numpy as np
    import h5py as h5
Unknown's avatar
Unknown committed
811
812

    main_h5_handle = h5.File(file_in_h5, 'r+')
Unknown's avatar
Unknown committed
813
    image_path = "/Frame_%04i/Channel_Finished" % img_num
Unknown's avatar
Unknown committed
814
    type_ref = type(main_h5_handle.get(image_path))
Unknown's avatar
Unknown committed
815
    temp = 1
Unknown's avatar
Unknown committed
816
    x = -1
Unknown's avatar
Unknown committed
817
818
819
    while temp:
        x = x + 1
        image_path = "%s/Filter_Step_%04i" % (image_path, x)
Unknown's avatar
Unknown committed
820
        image_temp = "%s/Lattice" % image_path
Unknown's avatar
Unknown committed
821
822
        temp_loc = main_h5_handle.get(image_temp)
        if type(temp_loc) == type_ref:
Unknown's avatar
Unknown committed
823
            temp = 0
Unknown's avatar
Unknown committed
824
825

    image_path = image_temp
Unknown's avatar
Unknown committed
826
    image_path = "%s/Positions" % image_path
Unknown's avatar
Unknown committed
827
828

    h5_image = main_h5_handle.get(image_path)
829
830
    pos = np.empty(h5_image.shape, dtype=h5_image.dtype)
    h5_image.read_direct(pos)
Unknown's avatar
Unknown committed
831

832
    return pos
Unknown's avatar
Unknown committed
833
834
835


def run_PCA_atoms(file_in_h5, img_num, box_width):
836
837
838
    import numpy as np
    import h5py as h5
    from scipy import linalg
Unknown's avatar
Unknown committed
839
840

    main_h5_handle = h5.File(file_in_h5, 'r+')
Unknown's avatar
Unknown committed
841
    image_path = "/Frame_%04i/Channel_Finished" % img_num
Unknown's avatar
Unknown committed
842
    type_ref = type(main_h5_handle.get(image_path))
Unknown's avatar
Unknown committed
843
    temp = 1
Unknown's avatar
Unknown committed
844
    x = -1
Unknown's avatar
Unknown committed
845
846
847
    while temp:
        x = x + 1
        image_path = "%s/Filter_Step_%04i" % (image_path, x)
Unknown's avatar
Unknown committed
848
        image_temp = "%s/Lattice" % image_path
Unknown's avatar
Unknown committed
849
850
        temp_loc = main_h5_handle.get(image_temp)
        if type(temp_loc) == type_ref:
Unknown's avatar
Unknown committed
851
            temp = 0
Unknown's avatar
Unknown committed
852
853

    pos_path = image_temp
Unknown's avatar
Unknown committed
854
    image_path = "%s/Positions" % pos_path
Unknown's avatar
Unknown committed
855
856

    h5_image = main_h5_handle.get(image_path)
857
858
    pos = np.empty(h5_image.shape, dtype=h5_image.dtype)
    h5_image.read_direct(pos)
Unknown's avatar
Unknown committed
859
860
861
    current_ref = h5_image.ref
    current_reg = h5_image.regionref[0:len(pos), 0:len(pos[0])]

Unknown's avatar
Unknown committed
862
    image_path = "/Frame_%04i/Channel_Current" % img_num
Unknown's avatar
Unknown committed
863
    image_path = "%s/Filter_Step_%04i" % (image_path, 0)
Unknown's avatar
Unknown committed
864
    image_path = "%s/Filtered_Image" % image_path
Unknown's avatar
Unknown committed
865
866

    h5_image = main_h5_handle.get(image_path)
867
868
    img2 = np.empty(h5_image.shape, dtype=h5_image.dtype)
    h5_image.read_direct(img2)
Unknown's avatar
Unknown committed
869

Unknown's avatar
Unknown committed
870
871
872
873
    sec_i_ref = h5_image.attrs.get("Spectroscopic_Indices")
    sec_i_reg = h5_image.attrs.get("Spectroscopic_Indices_Region")
    pos_i_ref = h5_image.attrs.get("Position_Indices")
    pos_i_reg = h5_image.attrs.get("Position_Indices_Region")
Unknown's avatar
Unknown committed
874

Unknown's avatar
Unknown committed
875
876
    spec_ind = main_h5_handle[sec_i_ref]
    spec_ind = spec_ind[sec_i_reg]
Unknown's avatar
Unknown committed
877

Unknown's avatar
Unknown committed
878
879
    posi_ind = main_h5_handle[pos_i_ref]
    posi_ind = posi_ind[pos_i_reg]
Unknown's avatar
Unknown committed
880
    img = np.empty([max(posi_ind[:, 0]), max(posi_ind[:, 1])], dtype=h5_image.dtype)
881
    try:
Unknown's avatar
Unknown committed
882
        img[posi_ind[:, 0] - 1, posi_ind[:, 1] - 1] = img2[0:len(img2), 0]
883
    except:
Unknown's avatar
Unknown committed
884
885
886
        img[posi_ind[:, 0] - 1, posi_ind[:, 1] - 1] = img2[0:len(img2)]

    sel_vec = np.zeros([len(pos[:, 0]), 1])
Unknown's avatar
Unknown committed
887
888
889
    kk = 0
    new_pos = []
    img_vectors = []
Unknown's avatar
Unknown committed
890
891
892
893
894
    for k1 in range(0, len(pos[:, 0])):
        if (box_width <= pos[k1, 0].round() <= len(img[:, 0]) - box_width and
                        box_width <= pos[k1, 1].round() <= len(img[0, :]) - box_width):
            sel_vec[k1] = 1
            new_pos.append(pos[k1, :])
Unknown's avatar
Unknown committed
895
            vector = img[pos[k1, 0] - box_width:pos[k1, 0] + box_width, pos[k1, 1] - box_width:pos[k1, 1] + box_width]
Unknown's avatar
Unknown committed
896
897
898
899
900
            img_vectors.append(vector.reshape([(box_width * 2) ** 2]))

    new_pos = array(new_pos)
    img_vectors = array(img_vectors)

901
    U, S, V = linalg.svd(img_vectors)
Unknown's avatar
Unknown committed
902
903
904
    V = V[0:len(S)]

    temp = 1
Unknown's avatar
Unknown committed
905
    x = -1
Unknown's avatar
Unknown committed
906
907
908
909
910
911

    while temp:
        x = x + 1
        image_temp = "%s/Analysis_%04i" % (pos_path, x)
        temp_loc = main_h5_handle.get(image_temp)
        if type(temp_loc) != type_ref:
Unknown's avatar
Unknown committed
912
            temp = 0
Unknown's avatar
Unknown committed
913

Unknown's avatar
Unknown committed
914
    image_path = image_temp
Unknown's avatar
Unknown committed
915
916
917
918
919
920
921

    [xx, yy] = np.meshgrid(range(1, len(U) + 1), range(1, len(U[0]) + 1))
    xx = np.reshape(xx, [1, len(U) * len(U[0])])
    yy = np.reshape(yy, [1, len(U) * len(U[0])])
    xy = [xx, yy]
    U = np.reshape(U, [1, len(U) * len(U[0])])

Unknown's avatar
Unknown committed
922
    image_path_sv = "%s/Spectroscopic_Values_01" % image_path
Unknown's avatar
Unknown committed
923
924
    main_h5_handle[image_path_sv] = [1]
    h5_image_new = main_h5_handle.get(image_path_sv)
Unknown's avatar
Unknown committed
925
926
927
    new_sv_ref = h5_image_new.ref
    new_sv_reg = h5_image_new.regionref[0]

Unknown's avatar
Unknown committed
928
    image_path_si = "%s/Spectroscopic_Indices_01" % image_path
Unknown's avatar
Unknown committed
929
930
    main_h5_handle[image_path_si] = [1]
    h5_image_new = main_h5_handle.get(image_path_si)
Unknown's avatar
Unknown committed
931
932
933
    new_si_ref = h5_image_new.ref
    new_si_reg = h5_image_new.regionref[0]

Unknown's avatar
Unknown committed
934
    image_path_sv = "%s/Position_Values_01" % image_path
Unknown's avatar
Unknown committed
935
936
    main_h5_handle[image_path_sv] = xy
    h5_image_new = main_h5_handle.get(image_path_sv)
Unknown's avatar
Unknown committed
937
938
939
    new_pv_ref = h5_image_new.ref
    new_pv_reg = h5_image_new.regionref[0:1, 0:len(xy[0])]

Unknown's avatar
Unknown committed
940
    image_path_si = "%s/Position_Indices_01" % image_path
Unknown's avatar
Unknown committed
941
942
    main_h5_handle[image_path_si] = xy
    h5_image_new = main_h5_handle.get(image_path_si)
Unknown's avatar
Unknown committed
943
944
945
    new_pi_ref = h5_image_new.ref
    new_pi_reg = h5_image_new.regionref[0:1, 0:len(xy[0])]

Unknown's avatar
Unknown committed
946
    image_path_b = "%s/Analysis_Data_01_U" % image_path
Unknown's avatar
Unknown committed
947
948
    main_h5_handle[image_path_b] = U
    h5_image_new = main_h5_handle.get(image_path_b)
949
    h5_new_attrs = h5_image_new.attrs
Unknown's avatar
Unknown committed
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
    h5_new_attrs["Filter_Name"] = "PCA_Atom_Shape"
    h5_new_attrs["Number_Of_Variables"] = 1
    h5_new_attrs["Variable_1_Name"] = "Box Width"
    h5_new_attrs["Variable_1_Value"] = box_width

    h5_new_attrs["Spectroscopic_Indices"] = new_si_ref
    h5_new_attrs["Spectroscopic_Indices_Region"] = new_si_reg
    h5_new_attrs["Spectroscopic_Values"] = new_sv_ref
    h5_new_attrs["Spectroscopic_Values_Region"] = new_sv_reg
    h5_new_attrs["Position_Indices"] = new_pi_ref
    h5_new_attrs["Position_Indices_Region"] = new_pi_reg
    h5_new_attrs["Position_Values"] = new_pv_ref
    h5_new_attrs["Position_Values_Region"] = new_pv_reg
    h5_new_attrs["Parent"] = current_ref
    h5_new_attrs["Parent_Region"] = current_reg
Unknown's avatar
Unknown committed
965
966
967

    [xy] = np.meshgrid(range(1, len(S) + 1))

Unknown's avatar
Unknown committed
968
    image_path_sv = "%s/Spectroscopic_Values_02" % image_path
Unknown's avatar
Unknown committed
969
970
    main_h5_handle[image_path_sv] = [1]
    h5_image_new = main_h5_handle.get(image_path_sv)
Unknown's avatar
Unknown committed
971
972
973
    new_sv_ref = h5_image_new.ref
    new_sv_reg = h5_image_new.regionref[0]

Unknown's avatar
Unknown committed
974
    image_path_si = "%s/Spectroscopic_Indices_02" % image_path
Unknown's avatar
Unknown committed
975
976
    main_h5_handle[image_path_si] = [1]
    h5_image_new = main_h5_handle.get(image_path_si)
Unknown's avatar
Unknown committed
977
978
979
    new_si_ref = h5_image_new.ref
    new_si_reg = h5_image_new.regionref[0]

Unknown's avatar
Unknown committed
980
    image_path_sv = "%s/Position_Values_02" % image_path
Unknown's avatar
Unknown committed
981
982
    main_h5_handle[image_path_sv] = xy
    h5_image_new = main_h5_handle.get(image_path_sv)
Unknown's avatar
Unknown committed
983
984
985
    new_pv_ref = h5_image_new.ref
    new_pv_reg = h5_image_new.regionref[0:len(xy)]