Mercurial > repos > rv43 > chess_tomo

comparison workflow/run_tomo.py @ 33:50923144fb56 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload for repository https://github.com/rolfverberg/galaxytools commit f8c4bdb31c20c468045ad5e6eb255a293244bc6c
author rv43
date Wed, 22 Mar 2023 16:46:26 +0000
parents 98eee43a8c8c
children
comparison
equal deleted inserted replaced
32:98eee43a8c8c 33:50923144fb56
13 except: 13 except:
14 pass 14 pass
15 15
16 from multiprocessing import cpu_count 16 from multiprocessing import cpu_count
17 from nexusformat.nexus import * 17 from nexusformat.nexus import *
18 from os import mkdir 18 from os import mkdir, environ
19 from os import path as os_path 19 from os import path as os_path
20 try: 20 try:
21 from skimage.transform import iradon 21 from skimage.transform import iradon
22 except: 22 except:
23 pass 23 pass
97 self.num_core = cpu_count() 97 self.num_core = cpu_count()
98 else: 98 else:
99 self.num_core = num_core 99 self.num_core = num_core
100 100
101 def __enter__(self): 101 def __enter__(self):
102 self.num_core_org = ne.set_num_threads(self.num_core) 102 self.num_core_org = ne.set_num_threads(min(self.num_core, ne.MAX_THREADS))
103 103
104 def __exit__(self, exc_type, exc_value, traceback): 104 def __exit__(self, exc_type, exc_value, traceback):
105 ne.set_num_threads(self.num_core_org) 105 ne.set_num_threads(self.num_core_org)
106 106
107 class Tomo: 107 class Tomo:
771 assert(len(tdf_stack) == 1) # TODO 771 assert(len(tdf_stack) == 1) # TODO
772 tdf_stack = tdf_stack[0] 772 tdf_stack = tdf_stack[0]
773 773
774 # Take median 774 # Take median
775 if tdf_stack.ndim == 2: 775 if tdf_stack.ndim == 2:
776 tdf = tdf_stack 776 tdf = tdf_stack.astype('float64')
777 elif tdf_stack.ndim == 3: 777 elif tdf_stack.ndim == 3:
778 tdf = np.median(tdf_stack, axis=0) 778 tdf = np.median(tdf_stack, axis=0)
779 del tdf_stack 779 del tdf_stack
780 else: 780 else:
781 raise ValueError(f'Invalid tdf_stack shape ({tdf_stack.shape})') 781 raise ValueError(f'Invalid tdf_stack shape ({tdf_stack.shape})')
839 corner of the frame where there is no sample but there is the direct X-ray 839 corner of the frame where there is no sample but there is the direct X-ray
840 beam because there is frame to frame fluctuations from the incoming beam. 840 beam because there is frame to frame fluctuations from the incoming beam.
841 We don’t typically account for them but potentially could. 841 We don’t typically account for them but potentially could.
842 """ 842 """
843 if tbf_stack.ndim == 2: 843 if tbf_stack.ndim == 2:
844 tbf = tbf_stack 844 tbf = tbf_stack.astype('float64')
845 elif tbf_stack.ndim == 3: 845 elif tbf_stack.ndim == 3:
846 tbf = np.median(tbf_stack, axis=0) 846 tbf = np.median(tbf_stack, axis=0)
847 del tbf_stack 847 del tbf_stack
848 else: 848 else:
849 raise ValueError(f'Invalid tbf_stack shape ({tbf_stacks.shape})') 849 raise ValueError(f'Invalid tbf_stack shape ({tbf_stacks.shape})')
854 else: 854 else:
855 logger.warning('Dark field unavailable') 855 logger.warning('Dark field unavailable')
856 856
857 # Set any non-positive values to one 857 # Set any non-positive values to one
858 # (avoid negative bright field values for spikes in dark field) 858 # (avoid negative bright field values for spikes in dark field)
859 tbf[tbf < 1] = 1 859 tbf[tbf < 1.0] = 1.0
860 860
861 # Plot bright field 861 # Plot bright field
862 if self.galaxy_flag: 862 if self.galaxy_flag:
863 quick_imshow(tbf, title='bright field', path='tomo_reduce_plots', 863 quick_imshow(tbf, title='bright field', path='tomo_reduce_plots',
864 save_fig=self.save_figs, save_only=self.save_only) 864 save_fig=self.save_figs, save_only=self.save_only)
1067 def _gen_tomo(self, nxentry, reduced_data): 1067 def _gen_tomo(self, nxentry, reduced_data):
1068 """Generate tomography fields. 1068 """Generate tomography fields.
1069 """ 1069 """
1070 # Get full bright field 1070 # Get full bright field
1071 tbf = np.asarray(reduced_data.data.bright_field) 1071 tbf = np.asarray(reduced_data.data.bright_field)
1072 tbf_shape = tbf.shape 1072 img_shape = tbf.shape
1073 1073
1074 # Get image bounds 1074 # Get image bounds
1075 img_x_bounds = tuple(reduced_data.get('img_x_bounds', (0, tbf_shape[0]))) 1075 img_x_bounds = tuple(reduced_data.get('img_x_bounds', (0, img_shape[0])))
1076 img_y_bounds = tuple(reduced_data.get('img_y_bounds', (0, tbf_shape[1]))) 1076 img_y_bounds = tuple(reduced_data.get('img_y_bounds', (0, img_shape[1])))
1077 if img_x_bounds == (0, img_shape[0]) and img_y_bounds == (0, img_shape[1]):
1078 resize_flag = False
1079 else:
1080 resize_flag = True
1077 1081
1078 # Get resized dark field 1082 # Get resized dark field
1079 # if 'dark_field' in data: 1083 if 'dark_field' in reduced_data.data:
1080 # tbf = np.asarray(reduced_data.data.dark_field[ 1084 if resize_flag:
1081 # img_x_bounds[0]:img_x_bounds[1],img_y_bounds[0]:img_y_bounds[1]]) 1085 tdf = np.asarray(reduced_data.data.dark_field[
1082 # else: 1086 img_x_bounds[0]:img_x_bounds[1],img_y_bounds[0]:img_y_bounds[1]])
1083 # logger.warning('Dark field unavailable') 1087 else:
1084 # tdf = None 1088 tdf = np.asarray(reduced_data.data.dark_field)
1085 tdf = None 1089 else:
1090 logger.warning('Dark field unavailable')
1091 tdf = None
1086 1092
1087 # Resize bright field 1093 # Resize bright field
1088 if img_x_bounds != (0, tbf.shape[0]) or img_y_bounds != (0, tbf.shape[1]): 1094 if resize_flag:
1089 tbf = tbf[img_x_bounds[0]:img_x_bounds[1],img_y_bounds[0]:img_y_bounds[1]] 1095 tbf = tbf[img_x_bounds[0]:img_x_bounds[1],img_y_bounds[0]:img_y_bounds[1]]
1090 1096
1091 # Get the tomography images 1097 # Get the tomography images
1092 image_key = nxentry.instrument.detector.get('image_key', None) 1098 image_key = nxentry.instrument.detector.get('image_key', None)
1093 if image_key and 'data' in nxentry.instrument.detector: 1099 if image_key and 'data' in nxentry.instrument.detector:
1143 if self.galaxy_flag: 1149 if self.galaxy_flag:
1144 path = 'tomo_reduce_plots' 1150 path = 'tomo_reduce_plots'
1145 else: 1151 else:
1146 path = self.output_folder 1152 path = self.output_folder
1147 for i, tomo_stack in enumerate(tomo_stacks): 1153 for i, tomo_stack in enumerate(tomo_stacks):
1148 # Resize the tomography images 1154 # Resize the tomography images as needed
1149 # Right now the range is the same for each set in the image stack. 1155 # Right now the range is the same for each set in the image stack
1150 if img_x_bounds != (0, tbf.shape[0]) or img_y_bounds != (0, tbf.shape[1]): 1156 if resize_flag:
1151 t0 = time() 1157 t0 = time()
1152 tomo_stack = tomo_stack[:,img_x_bounds[0]:img_x_bounds[1], 1158 tomo_stack = tomo_stack[:,img_x_bounds[0]:img_x_bounds[1],
1153 img_y_bounds[0]:img_y_bounds[1]].astype('float64') 1159 img_y_bounds[0]:img_y_bounds[1]].astype('float64')
1154 logger.debug(f'Resizing tomography images took {time()-t0:.2f} seconds') 1160 logger.debug(f'Resizing tomography images took {time()-t0:.2f} seconds')
1161 else:
1162 tomo_stack = tomo_stack.astype('float64')
1155 1163
1156 # Subtract dark field 1164 # Subtract dark field
1157 if tdf is not None: 1165 if tdf is not None:
1158 t0 = time() 1166 t0 = time()
1159 with set_numexpr_threads(self.num_core): 1167 with set_numexpr_threads(self.num_core):