Mercurial > repos > imgteam > projective_transformation

comparison projective_transformation.py @ 2:1ffdb07020ee 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/tools/projective_transformation/ commit 40b2eeecc9e15be0b710f5e5195ac6f3b3b83c39"
author imgteam
date Fri, 14 Jan 2022 03:07:08 +0000
parents 974cf4357707
children be9a815e2240
comparison
equal deleted inserted replaced
1:974cf4357707 2:1ffdb07020ee
1 import skimage.io 1 """
2 from skimage.transform import ProjectiveTransform 2 Copyright 2019-2022 Biomedical Computer Vision Group, Heidelberg University.
3 from scipy.ndimage import map_coordinates 3
4 Distributed under the MIT license.
5 See file LICENSE for detail or copy at https://opensource.org/licenses/MIT
6
7 """
8
9 import argparse
10 import imghdr
11 import warnings
12
4 import numpy as np 13 import numpy as np
5 import pandas as pd 14 import pandas as pd
6 import argparse 15 import skimage.color
7 import warnings 16 import skimage.io
8 import shutil 17 import tifffile
18 from scipy.ndimage import map_coordinates
19 from skimage.transform import ProjectiveTransform
9 20
10 21
11 def _stackcopy(a, b): 22 def _stackcopy(a, b):
12 if a.ndim == 3: 23 if a.ndim == 3:
13 a[:] = b[:, :, np.newaxis] 24 a[:] = b[:, :, np.newaxis]
22 coords_shape.append(shape[2]) 33 coords_shape.append(shape[2])
23 coords = np.empty(coords_shape, dtype=dtype) 34 coords = np.empty(coords_shape, dtype=dtype)
24 35
25 tf_coords = np.indices((cols, rows), dtype=dtype).reshape(2, -1).T 36 tf_coords = np.indices((cols, rows), dtype=dtype).reshape(2, -1).T
26 37
27 for i in range(0, (tf_coords.shape[0]//batch_size+1)): 38 for i in range(0, (tf_coords.shape[0] // batch_size + 1)):
28 tf_coords[batch_size*i:batch_size*(i+1)] = coord_map(tf_coords[batch_size*i:batch_size*(i+1)]) 39 tf_coords[batch_size * i:batch_size * (i + 1)] = coord_map(tf_coords[batch_size * i:batch_size * (i + 1)])
29 tf_coords = tf_coords.T.reshape((-1, cols, rows)).swapaxes(1, 2) 40 tf_coords = tf_coords.T.reshape((-1, cols, rows)).swapaxes(1, 2)
30 41
31 _stackcopy(coords[1, ...], tf_coords[0, ...]) 42 _stackcopy(coords[1, ...], tf_coords[0, ...])
32 _stackcopy(coords[0, ...], tf_coords[1, ...]) 43 _stackcopy(coords[0, ...], tf_coords[1, ...])
33 if len(shape) == 3: 44 if len(shape) == 3:
34 coords[2, ...] = range(shape[2]) 45 coords[2, ...] = range(shape[2])
35 46
36 return coords 47 return coords
37 48
38 49
39 def transform(moving_image, fixed_image, warp_matrix, out): 50 def transform(moving_fn, fixed_fn, warp_mat, output_fn):
40 moving_image = skimage.io.imread(moving_image)
41 fixed_image = skimage.io.imread(fixed_image)
42 warp_matrix = pd.read_csv(warp_matrix, delimiter="\t", header=None)
43 warp_matrix = np.array(warp_matrix)
44 51
45 trans = ProjectiveTransform(matrix=warp_matrix) 52 moving = skimage.io.imread(moving_fn)
46 warped_coords = warp_coords_batch(trans, fixed_image.shape) 53 extension = imghdr.what(moving_fn)
47 t = map_coordinates(moving_image, warped_coords, mode='reflect') 54 nDims = len(moving.shape)
55 assert nDims in [2, 3, 4, 5, 6], 'this tool only supports up to 6 dimensions'
56
57 if moving.shape[-1] in [3, 4] and nDims > 2:
58 isRGB = True
59 moving = np.transpose(moving, (nDims - 1,) + tuple(_ for _ in range(nDims - 1)))
60 else:
61 isRGB = False
62
63 if nDims > 3 or (nDims == 3 and not isRGB):
64 isMulCh = True
65 else:
66 isMulCh = False
67
68 fixed = skimage.io.imread(fixed_fn)
69 if fixed.shape[-1] in [3, 4] and len(fixed.shape) > 2:
70 hw_fixed = fixed.shape[-3:-1]
71 else:
72 hw_fixed = fixed.shape[-2:]
73
74 if isRGB or isMulCh:
75 shapeCh = moving.shape[0:-2]
76 nCh = np.prod(shapeCh)
77 moving = np.reshape(moving, (nCh,) + moving.shape[-2:])
78 warped_moving = np.zeros((nCh,) + hw_fixed, dtype=moving.dtype)
79
80 warp_mat = pd.read_csv(warp_mat, delimiter="\t", header=None)
81 warp_mat = np.array(warp_mat)
82 assert warp_mat.shape[0] in [3], 'only 2D image transformaton is supported'
83
84 trans = ProjectiveTransform(matrix=warp_mat)
85 warped_coords = warp_coords_batch(trans, hw_fixed)
86
87 if isMulCh or isRGB:
88 for i in range(nCh):
89 warped_moving[i, ...] = map_coordinates(moving[i, ...], warped_coords, cval=0)
90 warped_moving = np.reshape(warped_moving, shapeCh + warped_moving.shape[-2:])
91 if isRGB:
92 warped_moving = np.transpose(warped_moving, tuple(_ for _ in range(1, nDims)) + (0,))
93 else:
94 warped_moving = map_coordinates(moving, warped_coords, cval=0)
48 95
49 with warnings.catch_warnings(): 96 with warnings.catch_warnings():
50 warnings.simplefilter("ignore") 97 warnings.simplefilter("ignore")
51 skimage.io.imsave(out, t) 98 if isMulCh:
99 tifffile.imwrite(output_fn + '.tif', warped_moving, imagej=True, metadata={'mode': 'composite'})
100 else:
101 skimage.io.imsave(output_fn + '.' + extension, warped_moving)
52 102
53 103
54 if __name__ == "__main__": 104 if __name__ == "__main__":
55 parser = argparse.ArgumentParser(description="Transform the image") 105 parser = argparse.ArgumentParser(description="Transform the image")
56 parser.add_argument("fixed_image", help="Paste path to image.png that should be transformed") 106 parser.add_argument("fixed_image", help="Path to the fixed image")
57 parser.add_argument("moving_image", help="Paste path to fixed image.png") 107 parser.add_argument("moving_image", help="Path to the moving image (to be transformed)")
58 parser.add_argument("warp_matrix", help="Paste path to warp_matrix.csv that should be used for transformation") 108 parser.add_argument("warp_matrix", help="Path to the transformation matrix")
59 parser.add_argument("out", help="Paste path to file in which transformed image should be saved") 109 parser.add_argument("warped_image", help="Path to the output (transfirmed moving image)")
60 args = parser.parse_args() 110 args = parser.parse_args()
61 transform(args.moving_image, args.fixed_image, args.warp_matrix, args.out) 111 transform(args.moving_image, args.fixed_image, args.warp_matrix, args.warped_image)