comparison imagej2_analyze_skeleton_jython_script.py @ 1:1a7c29d5fc11 draft

"planemo upload for repository https://github.com/bgruening/galaxytools/tree/master/tools/image_processing/imagej2 commit 2afb24f3c81d625312186750a714d702363012b5"
author imgteam
date Mon, 28 Sep 2020 16:43:06 +0000
parents dc3935a0165a
children 666af1279007
comparison
equal deleted inserted replaced
0:dc3935a0165a 1:1a7c29d5fc11
1 import jython_utils
2 import math 1 import math
3 import sys 2 import sys
3
4 from ij import IJ 4 from ij import IJ
5 from sc.fiji.analyzeSkeleton import AnalyzeSkeleton_ 5 from sc.fiji.analyzeSkeleton import AnalyzeSkeleton_
6 6
7 BASIC_NAMES = [ 'Branches', 'Junctions', 'End-point Voxels', 7 BASIC_NAMES = ['Branches', 'Junctions', 'End-point Voxels',
8 'Junction Voxels', 'Slab Voxels', 'Average branch length', 8 'Junction Voxels', 'Slab Voxels', 'Average branch length',
9 'Triple Points', 'Quadruple Points', 'Maximum Branch Length' ] 9 'Triple Points', 'Quadruple Points', 'Maximum Branch Length']
10 DETAIL_NAMES = [ 'Skeleton ID', 'Branch length', 'V1 x', 'V1 y', 'V1 z', 'V2 x', 10 DETAIL_NAMES = ['Skeleton ID', 'Branch length', 'V1 x', 'V1 y', 'V1 z', 'V2 x',
11 'V2 y', 'V2 z', 'Euclidean distance' ] 11 'V2 y', 'V2 z', 'Euclidean distance']
12 OPTIONS = ['edm=Overwrite', 'iterations=1', 'count=1']
12 13
13 def get_euclidean_distance( vertex1, vertex2 ):
14 x1, y1, z1 = get_points( vertex1 )
15 x2, y2, z2 = get_points( vertex2 )
16 return math.sqrt( math.pow( ( x2 - x1 ), 2 ) +
17 math.pow( ( y2 - y1 ), 2 ) +
18 math.pow( ( z2 - z1 ), 2 ) )
19 14
20 def get_graph_length( graph ): 15 def get_euclidean_distance(vertex1, vertex2):
16 x1, y1, z1 = get_points(vertex1)
17 x2, y2, z2 = get_points(vertex2)
18 return math.sqrt(math.pow((x2 - x1), 2) + math.pow((y2 - y1), 2) + math.pow((z2 - z1), 2))
19
20
21 def get_graph_length(graph):
21 length = 0 22 length = 0
22 for edge in graph.getEdges(): 23 for edge in graph.getEdges():
23 length = length + edge.getLength() 24 length = length + edge.getLength()
24 return length 25 return length
25 26
26 def get_points( vertex ): 27
28 def get_points(vertex):
27 # An array of Point, which has attributes x,y,z. 29 # An array of Point, which has attributes x,y,z.
28 point = vertex.getPoints()[ 0 ] 30 point = vertex.getPoints()[0]
29 return point.x, point.y, point.z 31 return point.x, point.y, point.z
30 32
31 def get_sorted_edge_lengths( graph ): 33
34 def get_sorted_edge_lengths(graph):
32 # Return graph edges sorted from longest to shortest. 35 # Return graph edges sorted from longest to shortest.
33 edges = graph.getEdges() 36 edges = graph.getEdges()
34 edges = sorted( edges, key=lambda edge: edge.getLength(), reverse=True ) 37 edges = sorted(edges, key=lambda edge: edge.getLength(), reverse=True)
35 return edges 38 return edges
36 39
37 def get_sorted_graph_lengths( result ): 40
41 def get_sorted_graph_lengths(result):
38 # Get the separate graphs (skeletons). 42 # Get the separate graphs (skeletons).
39 graphs = result.getGraph() 43 graphs = result.getGraph()
40 # Sort graphs from longest to shortest. 44 # Sort graphs from longest to shortest.
41 graphs = sorted( graphs, key=lambda g: get_graph_length( g ), reverse=True ) 45 graphs = sorted(graphs, key=lambda g: get_graph_length(g), reverse=True)
42 return graphs 46 return graphs
43 47
44 def save( result, output, show_detailed_info, calculate_largest_shortest_path, sep='\t' ): 48
45 num_trees = int( result.getNumOfTrees() ) 49 def save(result, output, show_detailed_info, calculate_largest_shortest_path, sep='\t'):
46 outf = open( output, 'wb' ) 50 num_trees = int(result.getNumOfTrees())
47 outf.write( '# %s\n' % sep.join( BASIC_NAMES ) ) 51 outf = open(output, 'wb')
48 for index in range( num_trees ): 52 outf.write('# %s\n' % sep.join(BASIC_NAMES))
49 outf.write( '%d%s' % ( result.getBranches()[ index ], sep ) ) 53 for index in range(num_trees):
50 outf.write( '%d%s' % ( result.getJunctions()[ index ], sep ) ) 54 outf.write('%d%s' % (result.getBranches()[index], sep))
51 outf.write( '%d%s' % ( result.getEndPoints()[ index ], sep ) ) 55 outf.write('%d%s' % (result.getJunctions()[index], sep))
52 outf.write( '%d%s' % ( result.getJunctionVoxels()[ index ], sep ) ) 56 outf.write('%d%s' % (result.getEndPoints()[index], sep))
53 outf.write( '%d%s' % ( result.getSlabs()[ index ], sep ) ) 57 outf.write('%d%s' % (result.getJunctionVoxels()[index], sep))
54 outf.write( '%.3f%s' % ( result.getAverageBranchLength()[ index ], sep ) ) 58 outf.write('%d%s' % (result.getSlabs()[index], sep))
55 outf.write( '%d%s' % ( result.getTriples()[ index ], sep ) ) 59 outf.write('%.3f%s' % (result.getAverageBranchLength()[index], sep))
56 outf.write( '%d%s' % ( result.getQuadruples()[ index ], sep ) ) 60 outf.write('%d%s' % (result.getTriples()[index], sep))
57 outf.write( '%.3f' % result.getMaximumBranchLength()[ index ] ) 61 outf.write('%d%s' % (result.getQuadruples()[index], sep))
62 outf.write('%.3f' % result.getMaximumBranchLength()[index])
58 if calculate_largest_shortest_path: 63 if calculate_largest_shortest_path:
59 outf.write( '%s%.3f%s' % ( sep, result.shortestPathList.get( index ), sep ) ) 64 outf.write('%s%.3f%s' % (sep, result.shortestPathList.get(index), sep))
60 outf.write( '%d%s' % ( result.spStartPosition[ index ][ 0 ], sep ) ) 65 outf.write('%d%s' % (result.spStartPosition[index][0], sep))
61 outf.write( '%d%s' % ( result.spStartPosition[ index ][ 1 ], sep ) ) 66 outf.write('%d%s' % (result.spStartPosition[index][1], sep))
62 outf.write( '%d\n' % result.spStartPosition[ index ][ 2 ] ) 67 outf.write('%d\n' % result.spStartPosition[index][2])
63 else: 68 else:
64 outf.write( '\n' ) 69 outf.write('\n')
65 if show_detailed_info: 70 if show_detailed_info:
66 outf.write( '# %s\n' % sep.join( DETAIL_NAMES ) ) 71 outf.write('# %s\n' % sep.join(DETAIL_NAMES))
67 # The following index is a placeholder for the skeleton ID. 72 # The following index is a placeholder for the skeleton ID.
68 # The terms "graph" and "skeleton" refer to the same thing. 73 # The terms "graph" and "skeleton" refer to the same thing.
69 # Also, the SkeletonResult.java code states that the 74 # Also, the SkeletonResult.java code states that the
70 # private Graph[] graph object is an array of graphs (one 75 # private Graph[] graph object is an array of graphs (one
71 # per tree). 76 # per tree).
72 for index, graph in enumerate( get_sorted_graph_lengths( result ) ): 77 for index, graph in enumerate(get_sorted_graph_lengths(result)):
73 for edge in get_sorted_edge_lengths( graph ): 78 for edge in get_sorted_edge_lengths(graph):
74 vertex1 = edge.getV1() 79 vertex1 = edge.getV1()
75 x1, y1, z1 = get_points( vertex1 ) 80 x1, y1, z1 = get_points(vertex1)
76 vertex2 = edge.getV2() 81 vertex2 = edge.getV2()
77 x2, y2, z2 = get_points( vertex2 ) 82 x2, y2, z2 = get_points(vertex2)
78 outf.write( '%d%s' % ( index+1, sep ) ) 83 outf.write('%d%s' % (index + 1, sep))
79 outf.write( '%.3f%s' % ( edge.getLength(), sep ) ) 84 outf.write('%.3f%s' % (edge.getLength(), sep))
80 outf.write( '%d%s' % ( x1, sep ) ) 85 outf.write('%d%s' % (x1, sep))
81 outf.write( '%d%s' % ( y1, sep ) ) 86 outf.write('%d%s' % (y1, sep))
82 outf.write( '%d%s' % ( z1, sep ) ) 87 outf.write('%d%s' % (z1, sep))
83 outf.write( '%d%s' % ( x2, sep ) ) 88 outf.write('%d%s' % (x2, sep))
84 outf.write( '%d%s' % ( y2, sep ) ) 89 outf.write('%d%s' % (y2, sep))
85 outf.write( '%d%s' % ( z2, sep ) ) 90 outf.write('%d%s' % (z2, sep))
86 outf.write( '%.3f' % get_euclidean_distance( vertex1, vertex2 ) ) 91 outf.write('%.3f' % get_euclidean_distance(vertex1, vertex2))
87 if calculate_largest_shortest_path: 92 if calculate_largest_shortest_path:
88 # Keep number of separated items the same for each line. 93 # Keep number of separated items the same for each line.
89 outf.write( '%s %s' % ( sep, sep ) ) 94 outf.write('%s %s' % (sep, sep))
90 outf.write( ' %s' % sep ) 95 outf.write(' %s' % sep)
91 outf.write( ' %s' % sep ) 96 outf.write(' %s' % sep)
92 outf.write( ' \n' ) 97 outf.write(' \n')
93 else: 98 else:
94 outf.write( '\n' ) 99 outf.write('\n')
95 outf.close() 100 outf.close()
101
96 102
97 # Fiji Jython interpreter implements Python 2.5 which does not 103 # Fiji Jython interpreter implements Python 2.5 which does not
98 # provide support for argparse. 104 # provide support for argparse.
99 error_log = sys.argv[ -8 ] 105 error_log = sys.argv[-8]
100 input = sys.argv[ -7 ] 106 input = sys.argv[-7]
101 black_background = jython_utils.asbool( sys.argv[ -6 ] ) 107 black_background = sys.argv[-6] == "yes"
102 prune_cycle_method = sys.argv[ -5 ] 108 prune_cycle_method = sys.argv[-5]
103 prune_ends = jython_utils.asbool( sys.argv[ -4 ] ) 109 prune_ends = sys.argv[-4] == "yes"
104 calculate_largest_shortest_path = jython_utils.asbool( sys.argv[ -3 ] ) 110 calculate_largest_shortest_path = sys.argv[-3] == "yes"
105 if calculate_largest_shortest_path: 111 if calculate_largest_shortest_path:
106 BASIC_NAMES.extend( [ 'Longest Shortest Path', 'spx', 'spy', 'spz' ] ) 112 BASIC_NAMES.extend(['Longest Shortest Path', 'spx', 'spy', 'spz'])
107 DETAIL_NAMES.extend( [ ' ', ' ', ' ', ' ' ] ) 113 DETAIL_NAMES.extend([' ', ' ', ' ', ' '])
108 show_detailed_info = jython_utils.asbool( sys.argv[ -2 ] ) 114 show_detailed_info = sys.argv[-2] == "yes"
109 output = sys.argv[ -1 ] 115 output = sys.argv[-1]
110 116
111 # Open the input image file. 117 # Open the input image file.
112 input_image_plus = IJ.openImage( input ) 118 input_image_plus = IJ.openImage(input)
113 119
114 # Create a copy of the image. 120 # Create a copy of the image.
115 input_image_plus_copy = input_image_plus.duplicate() 121 input_image_plus_copy = input_image_plus.duplicate()
116 image_processor_copy = input_image_plus_copy.getProcessor() 122 image_processor_copy = input_image_plus_copy.getProcessor()
117 123
118 try: 124 # Set binary options.
119 # Set binary options. 125 options_list = OPTIONS
120 options = jython_utils.get_binary_options( black_background=black_background ) 126 if black_background:
121 IJ.run( input_image_plus_copy, "Options...", options ) 127 options_list.append("black")
128 options = " ".join(options_list)
129 IJ.run(input_image_plus_copy, "Options...", options)
122 130
123 # Convert image to binary if necessary. 131 # Convert image to binary if necessary.
124 if not image_processor_copy.isBinary(): 132 if not image_processor_copy.isBinary():
125 IJ.run( input_image_plus_copy, "Make Binary", "" ) 133 IJ.run(input_image_plus_copy, "Make Binary", "")
126 134
127 # Run AnalyzeSkeleton 135 # Run AnalyzeSkeleton
128 analyze_skeleton = AnalyzeSkeleton_() 136 analyze_skeleton = AnalyzeSkeleton_()
129 analyze_skeleton.setup( "", input_image_plus_copy ) 137 analyze_skeleton.setup("", input_image_plus_copy)
130 if prune_cycle_method == 'none': 138 if prune_cycle_method == 'none':
131 prune_index = analyze_skeleton.NONE 139 prune_index = analyze_skeleton.NONE
132 elif prune_cycle_method == 'shortest_branch': 140 elif prune_cycle_method == 'shortest_branch':
133 prune_index = analyze_skeleton.SHORTEST_BRANCH 141 prune_index = analyze_skeleton.SHORTEST_BRANCH
134 elif prune_cycle_method == 'lowest_intensity_voxel': 142 elif prune_cycle_method == 'lowest_intensity_voxel':
135 prune_index = analyze_skeleton.LOWEST_INTENSITY_VOXEL 143 prune_index = analyze_skeleton.LOWEST_INTENSITY_VOXEL
136 elif prune_cycle_method == 'lowest_intensity_branch': 144 elif prune_cycle_method == 'lowest_intensity_branch':
137 prune_index = analyze_skeleton.LOWEST_INTENSITY_BRANCH 145 prune_index = analyze_skeleton.LOWEST_INTENSITY_BRANCH
138 result = analyze_skeleton.run( prune_index, 146 result = analyze_skeleton.run(prune_index, prune_ends, calculate_largest_shortest_path, input_image_plus_copy, True, True)
139 prune_ends, 147 # Save the results.
140 calculate_largest_shortest_path, 148 save(result, output, show_detailed_info, calculate_largest_shortest_path)
141 input_image_plus_copy,
142 True,
143 True )
144 # Save the results.
145 save( result, output, show_detailed_info, calculate_largest_shortest_path )
146 except Exception, e:
147 jython_utils.handle_error( error_log, str( e ) )