Mercurial > repos > imgteam > spot_detection_2d

comparison spot_detection_2d.py @ 1:859dd1c11ac0 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
"planemo upload for repository https://github.com/BMCV/galaxy-image-analysis/tree/master/tools/spot_detection_2d/ commit 9f103372d66ae7e3c5c385bd444b2a80e51cdae6"
author imgteam
date Sat, 19 Feb 2022 20:43:29 +0000
parents d78372040976
children 4645b356bf3b
comparison
equal deleted inserted replaced
0:d78372040976 1:859dd1c11ac0
11 11
12 import imageio 12 import imageio
13 import numpy as np 13 import numpy as np
14 import pandas as pd 14 import pandas as pd
15 from skimage.feature import peak_local_max 15 from skimage.feature import peak_local_max
16 from skimage.filters import gaussian 16 from skimage.filters import gaussian, laplace
17 17
18 18
19 def getbr(xy, img, nb, firstn): 19 def getbr(xy, img, nb, firstn):
20 ndata = xy.shape[0] 20 ndata = xy.shape[0]
21 br = np.empty((ndata, 1)) 21 br = np.empty((ndata, 1))
25 timg = img[xy[j, 1] - nb - 1:xy[j, 1] + nb, xy[j, 0] - nb - 1:xy[j, 0] + nb] 25 timg = img[xy[j, 1] - nb - 1:xy[j, 1] + nb, xy[j, 0] - nb - 1:xy[j, 0] + nb]
26 br[j] = np.mean(np.sort(timg, axis=None)[-firstn:]) 26 br[j] = np.mean(np.sort(timg, axis=None)[-firstn:])
27 return br 27 return br
28 28
29 29
30 def spot_detection(fn_in, fn_out, frame_1st=1, frame_end=0, typ_br='smoothed', th=10, ssig=1, bd=10): 30 def spot_detection(fn_in, fn_out, frame_1st=1, frame_end=0,
31 typ_filter='Gauss', ssig=1, th=10,
32 typ_br='smoothed', bd=10):
31 ims_ori = imageio.mimread(fn_in, format='TIFF') 33 ims_ori = imageio.mimread(fn_in, format='TIFF')
32 ims_smd = np.zeros((len(ims_ori), ims_ori[0].shape[0], ims_ori[0].shape[1]), dtype='float64') 34 ims_smd = np.zeros((len(ims_ori), ims_ori[0].shape[0], ims_ori[0].shape[1]), dtype='float64')
33 if frame_end == 0 or frame_end > len(ims_ori): 35 if frame_end == 0 or frame_end > len(ims_ori):
34 frame_end = len(ims_ori) 36 frame_end = len(ims_ori)
35 37
38 ims_smd_max = np.max(ims_smd) 40 ims_smd_max = np.max(ims_smd)
39 41
40 txyb_all = np.array([]).reshape(0, 4) 42 txyb_all = np.array([]).reshape(0, 4)
41 for i in range(frame_1st - 1, frame_end): 43 for i in range(frame_1st - 1, frame_end):
42 tmp = np.copy(ims_smd[i, :, :]) 44 tmp = np.copy(ims_smd[i, :, :])
45 if typ_filter == 'LoG':
46 tmp = laplace(tmp)
47
43 tmp[tmp < th * ims_smd_max / 100] = 0 48 tmp[tmp < th * ims_smd_max / 100] = 0
44 coords = peak_local_max(tmp, min_distance=1) 49 coords = peak_local_max(tmp, min_distance=1)
45 idx_to_del = np.where((coords[:, 0] <= bd) | (coords[:, 0] >= tmp.shape[0] - bd) | 50 idx_to_del = np.where((coords[:, 0] <= bd) | (coords[:, 0] >= tmp.shape[0] - bd) |
46 (coords[:, 1] <= bd) | (coords[:, 1] >= tmp.shape[1] - bd)) 51 (coords[:, 1] <= bd) | (coords[:, 1] >= tmp.shape[1] - bd))
47 coords = np.delete(coords, idx_to_del[0], axis=0) 52 coords = np.delete(coords, idx_to_del[0], axis=0)
64 df['INTENSITY'] = txyb_all[:, 3] 69 df['INTENSITY'] = txyb_all[:, 3]
65 df.to_csv(fn_out, index=False, float_format='%.2f', sep="\t") 70 df.to_csv(fn_out, index=False, float_format='%.2f', sep="\t")
66 71
67 72
68 if __name__ == "__main__": 73 if __name__ == "__main__":
69 parser = argparse.ArgumentParser(description="Spot detection based on local maxima") 74 parser = argparse.ArgumentParser(description="Spot detection")
70 parser.add_argument("fn_in", help="Name of input image sequence (stack)") 75 parser.add_argument("fn_in", help="Name of input image sequence (stack)")
71 parser.add_argument("fn_out", help="Name of output file to save the coordinates and intensities of detected spots") 76 parser.add_argument("fn_out", help="Name of output file to save the coordinates and intensities of detected spots")
72 parser.add_argument("frame_1st", type=int, help="Index for the starting frame to detect spots (1 for first frame of the stack)") 77 parser.add_argument("frame_1st", type=int, help="Index for the starting frame to detect spots (1 for first frame of the stack)")
73 parser.add_argument("frame_end", type=int, help="Index for the last frame to detect spots (0 for the last frame of the stack)") 78 parser.add_argument("frame_end", type=int, help="Index for the last frame to detect spots (0 for the last frame of the stack)")
79 parser.add_argument("filter", help="Detection filter")
80 parser.add_argument("ssig", type=float, help="Sigma of the Gaussian for noise suppression")
81 parser.add_argument("thres", type=float, help="Percentage of the global maximal for thresholding candidate spots")
74 parser.add_argument("typ_intens", help="smoothed or robust (for measuring the intensities of spots)") 82 parser.add_argument("typ_intens", help="smoothed or robust (for measuring the intensities of spots)")
75 parser.add_argument("thres", type=float, help="Percentage of the global maximal intensity for thresholding candidate spots")
76 parser.add_argument("ssig", type=float, help="Sigma of the Gaussian filter for noise suppression")
77 parser.add_argument("bndy", type=int, help="Number of pixels (Spots close to image boundaries will be ignored)") 83 parser.add_argument("bndy", type=int, help="Number of pixels (Spots close to image boundaries will be ignored)")
78 args = parser.parse_args() 84 args = parser.parse_args()
79 spot_detection(args.fn_in, 85 spot_detection(args.fn_in, args.fn_out,
80 args.fn_out, 86 frame_1st=args.frame_1st, frame_end=args.frame_end,
81 frame_1st=args.frame_1st, 87 typ_filter=args.filter, ssig=args.ssig, th=args.thres,
82 frame_end=args.frame_end, 88 typ_br=args.typ_intens, bd=args.bndy)
83 typ_br=args.typ_intens,
84 th=args.thres,
85 ssig=args.ssig,
86 bd=args.bndy)