comparison workflow/run_tomo.py @ 75:d5e1d4ea2b7e draft default tip

planemo upload for repository https://github.com/rolfverberg/galaxytools commit 6afde341a94586fe3972bdbbfbf5dabd5e8dec69
author rv43
date Thu, 23 Mar 2023 13:39:14 +0000
parents 1cf15b61cd83
children
comparison
equal deleted inserted replaced
74:4f4ee8db5f67 75:d5e1d4ea2b7e
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:
162 logger.warning(f'num_core = {self.num_core} is larger than the number of available ' 162 logger.warning(f'num_core = {self.num_core} is larger than the number of available '
163 f'processors and reduced to {cpu_count()}') 163 f'processors and reduced to {cpu_count()}')
164 self.num_core= cpu_count() 164 self.num_core= cpu_count()
165 165
166 def read(self, filename): 166 def read(self, filename):
167 extension = os_path.splitext(filename)[1] 167 logger.info(f'looking for {filename}')
168 if extension == '.yml' or extension == '.yaml': 168 if self.galaxy_flag:
169 with open(filename, 'r') as f: 169 try:
170 config = safe_load(f) 170 with open(filename, 'r') as f:
171 # if len(config) > 1: 171 config = safe_load(f)
172 # raise ValueError(f'Multiple root entries in {filename} not yet implemented') 172 return(config)
173 # if len(list(config.values())[0]) > 1: 173 except:
174 # raise ValueError(f'Multiple sample maps in {filename} not yet implemented') 174 try:
175 return(config) 175 with NXFile(filename, mode='r') as nxfile:
176 elif extension == '.nxs': 176 nxroot = nxfile.readfile()
177 with NXFile(filename, mode='r') as nxfile: 177 return(nxroot)
178 nxroot = nxfile.readfile() 178 except:
179 return(nxroot) 179 raise ValueError(f'Unable to open ({filename})')
180 else: 180 else:
181 raise ValueError(f'Invalid filename extension ({extension})') 181 extension = os_path.splitext(filename)[1]
182 if extension == '.yml' or extension == '.yaml':
183 with open(filename, 'r') as f:
184 config = safe_load(f)
185 # if len(config) > 1:
186 # raise ValueError(f'Multiple root entries in {filename} not yet implemented')
187 # if len(list(config.values())[0]) > 1:
188 # raise ValueError(f'Multiple sample maps in {filename} not yet implemented')
189 return(config)
190 elif extension == '.nxs':
191 with NXFile(filename, mode='r') as nxfile:
192 nxroot = nxfile.readfile()
193 return(nxroot)
194 else:
195 raise ValueError(f'Invalid filename extension ({extension})')
182 196
183 def write(self, data, filename): 197 def write(self, data, filename):
184 extension = os_path.splitext(filename)[1] 198 extension = os_path.splitext(filename)[1]
185 if extension == '.yml' or extension == '.yaml': 199 if extension == '.yml' or extension == '.yaml':
186 with open(filename, 'w') as f: 200 with open(filename, 'w') as f:
187 safe_dump(data, f) 201 safe_dump(data, f)
188 elif extension == '.nxs': 202 elif extension == '.nxs' or extension == '.nex':
189 data.save(filename, mode='w') 203 data.save(filename, mode='w')
190 elif extension == '.nc': 204 elif extension == '.nc':
191 data.to_netcdf(os_path=filename) 205 data.to_netcdf(os_path=filename)
192 else: 206 else:
193 raise ValueError(f'Invalid filename extension ({extension})') 207 raise ValueError(f'Invalid filename extension ({extension})')
285 if not isinstance(nxroot, NXroot): 299 if not isinstance(nxroot, NXroot):
286 raise ValueError(f'Invalid parameter nxroot ({nxroot})') 300 raise ValueError(f'Invalid parameter nxroot ({nxroot})')
287 nxentry = nxroot[nxroot.attrs['default']] 301 nxentry = nxroot[nxroot.attrs['default']]
288 if not isinstance(nxentry, NXentry): 302 if not isinstance(nxentry, NXentry):
289 raise ValueError(f'Invalid nxentry ({nxentry})') 303 raise ValueError(f'Invalid nxentry ({nxentry})')
290 if self.galaxy_flag: 304 if center_rows is not None:
291 if center_rows is not None: 305 if self.galaxy_flag:
292 center_rows = tuple(center_rows) 306 if not is_int_pair(center_rows, ge=-1):
293 if not is_int_pair(center_rows):
294 raise ValueError(f'Invalid parameter center_rows ({center_rows})') 307 raise ValueError(f'Invalid parameter center_rows ({center_rows})')
295 elif center_rows is not None: 308 if (center_rows[0] != -1 and center_rows[1] != -1 and
296 logger.warning(f'Ignoring parameter center_rows ({center_rows})') 309 center_rows[0] > center_rows[1]):
297 center_rows = None 310 center_rows = (center_rows[1], center_rows[0])
311 else:
312 center_rows = tuple(center_rows)
313 else:
314 logger.warning(f'Ignoring parameter center_rows ({center_rows})')
315 center_rows = None
298 if self.galaxy_flag: 316 if self.galaxy_flag:
299 if center_stack_index is not None and not is_int(center_stack_index, ge=0): 317 if center_stack_index is not None and not is_int(center_stack_index, ge=0):
300 raise ValueError(f'Invalid parameter center_stack_index ({center_stack_index})') 318 raise ValueError(f'Invalid parameter center_stack_index ({center_stack_index})')
301 elif center_stack_index is not None: 319 elif center_stack_index is not None:
302 logger.warning(f'Ignoring parameter center_stack_index ({center_stack_index})') 320 logger.warning(f'Ignoring parameter center_stack_index ({center_stack_index})')
358 376
359 # Lower row center 377 # Lower row center
360 if self.test_mode: 378 if self.test_mode:
361 lower_row = self.test_config['lower_row'] 379 lower_row = self.test_config['lower_row']
362 elif self.galaxy_flag: 380 elif self.galaxy_flag:
363 if center_rows is None: 381 if center_rows is None or center_rows[0] == -1:
364 lower_row = 0 382 lower_row = 0
365 else: 383 else:
366 lower_row = min(center_rows) 384 lower_row = center_rows[0]
367 if not 0 <= lower_row < tomo_fields_shape[2]-1: 385 if not 0 <= lower_row < tomo_fields_shape[2]-1:
368 raise ValueError(f'Invalid parameter center_rows ({center_rows})') 386 raise ValueError(f'Invalid parameter center_rows ({center_rows})')
369 else: 387 else:
370 lower_row = select_one_image_bound( 388 lower_row = select_one_image_bound(
371 nxentry.reduced_data.data.tomo_fields[center_stack_index,0,:,:], 0, bound=0, 389 nxentry.reduced_data.data.tomo_fields[center_stack_index,0,:,:], 0, bound=0,
372 title=f'theta={round(thetas[0], 2)+0}', 390 title=f'theta={round(thetas[0], 2)+0}',
373 bound_name='row index to find lower center', default=default, raise_error=True) 391 bound_name='row index to find lower center', default=default, raise_error=True)
374 logger.debug('Finding center...') 392 logger.debug('Finding center...')
375 t0 = time() 393 t0 = time()
376 lower_center_offset = self._find_center_one_plane( 394 lower_center_offset = self._find_center_one_plane(
377 #np.asarray(nxentry.reduced_data.data.tomo_fields[center_stack_index,:,lower_row,:]),
378 nxentry.reduced_data.data.tomo_fields[center_stack_index,:,lower_row,:], 395 nxentry.reduced_data.data.tomo_fields[center_stack_index,:,lower_row,:],
379 lower_row, thetas, eff_pixel_size, cross_sectional_dim, path=path, 396 lower_row, thetas, eff_pixel_size, cross_sectional_dim, path=path,
380 num_core=self.num_core) 397 num_core=self.num_core)
381 logger.debug(f'... done in {time()-t0:.2f} seconds') 398 logger.debug(f'... done in {time()-t0:.2f} seconds')
382 logger.debug(f'lower_row = {lower_row:.2f}') 399 logger.debug(f'lower_row = {lower_row:.2f}')
384 401
385 # Upper row center 402 # Upper row center
386 if self.test_mode: 403 if self.test_mode:
387 upper_row = self.test_config['upper_row'] 404 upper_row = self.test_config['upper_row']
388 elif self.galaxy_flag: 405 elif self.galaxy_flag:
389 if center_rows is None: 406 if center_rows is None or center_rows[1] == -1:
390 upper_row = tomo_fields_shape[2]-1 407 upper_row = tomo_fields_shape[2]-1
391 else: 408 else:
392 upper_row = max(center_rows) 409 upper_row = center_rows[1]
393 if not lower_row < upper_row < tomo_fields_shape[2]: 410 if not lower_row < upper_row < tomo_fields_shape[2]:
394 raise ValueError(f'Invalid parameter center_rows ({center_rows})') 411 raise ValueError(f'Invalid parameter center_rows ({center_rows})')
395 else: 412 else:
396 upper_row = select_one_image_bound( 413 upper_row = select_one_image_bound(
397 nxentry.reduced_data.data.tomo_fields[center_stack_index,0,:,:], 0, 414 nxentry.reduced_data.data.tomo_fields[center_stack_index,0,:,:], 0,
497 tomo_recon_stacks[i] = tomo_recon_stack 514 tomo_recon_stacks[i] = tomo_recon_stack
498 515
499 # Resize the reconstructed tomography data 516 # Resize the reconstructed tomography data
500 # reconstructed data order in each stack: row/z,x,y 517 # reconstructed data order in each stack: row/z,x,y
501 if self.test_mode: 518 if self.test_mode:
502 x_bounds = self.test_config.get('x_bounds') 519 x_bounds = tuple(self.test_config.get('x_bounds'))
503 y_bounds = self.test_config.get('y_bounds') 520 y_bounds = tuple(self.test_config.get('y_bounds'))
504 z_bounds = None 521 z_bounds = None
505 elif self.galaxy_flag: 522 elif self.galaxy_flag:
506 if x_bounds is not None and not is_int_pair(x_bounds, ge=0, 523 x_max = tomo_recon_stacks[0].shape[1]
507 lt=tomo_recon_stacks[0].shape[1]): 524 if x_bounds is None:
525 x_bounds = (0, x_max)
526 elif is_int_pair(x_bounds, ge=-1, le=x_max):
527 x_bounds = tuple(x_bounds)
528 if x_bounds[0] == -1:
529 x_bounds = (0, x_bounds[1])
530 if x_bounds[1] == -1:
531 x_bounds = (x_bounds[0], x_max)
532 if not is_index_range(x_bounds, ge=0, le=x_max):
508 raise ValueError(f'Invalid parameter x_bounds ({x_bounds})') 533 raise ValueError(f'Invalid parameter x_bounds ({x_bounds})')
509 if y_bounds is not None and not is_int_pair(y_bounds, ge=0, 534 y_max = tomo_recon_stacks[0].shape[1]
510 lt=tomo_recon_stacks[0].shape[1]): 535 if y_bounds is None:
536 y_bounds = (0, y_max)
537 elif is_int_pair(y_bounds, ge=-1, le=y_max):
538 y_bounds = tuple(y_bounds)
539 if y_bounds[0] == -1:
540 y_bounds = (0, y_bounds[1])
541 if y_bounds[1] == -1:
542 y_bounds = (y_bounds[0], y_max)
543 if not is_index_range(y_bounds, ge=0, le=y_max):
511 raise ValueError(f'Invalid parameter y_bounds ({y_bounds})') 544 raise ValueError(f'Invalid parameter y_bounds ({y_bounds})')
512 z_bounds = None 545 z_bounds = None
513 else: 546 else:
514 x_bounds, y_bounds, z_bounds = self._resize_reconstructed_data(tomo_recon_stacks) 547 x_bounds, y_bounds, z_bounds = self._resize_reconstructed_data(tomo_recon_stacks)
515 if x_bounds is None: 548 if x_bounds is None:
533 566
534 # Plot a few reconstructed image slices 567 # Plot a few reconstructed image slices
535 if num_tomo_stacks == 1: 568 if num_tomo_stacks == 1:
536 basetitle = 'recon' 569 basetitle = 'recon'
537 else: 570 else:
538 basetitle = f'recon stack {i+1}' 571 basetitle = f'recon stack'
539 for i, stack in enumerate(tomo_recon_stacks): 572 for i, stack in enumerate(tomo_recon_stacks):
540 title = f'{basetitle} {res_title} xslice{x_slice}' 573 title = f'{basetitle} {i+1} {res_title} xslice{x_slice}'
541 quick_imshow(stack[z_range[0]:z_range[1],x_slice,y_range[0]:y_range[1]], 574 quick_imshow(stack[z_range[0]:z_range[1],x_slice,y_range[0]:y_range[1]],
542 title=title, path=path, save_fig=self.save_figs, save_only=self.save_only, 575 title=title, path=path, save_fig=self.save_figs, save_only=self.save_only,
543 block=self.block) 576 block=self.block)
544 title = f'{basetitle} {res_title} yslice{y_slice}' 577 title = f'{basetitle} {i+1} {res_title} yslice{y_slice}'
545 quick_imshow(stack[z_range[0]:z_range[1],x_range[0]:x_range[1],y_slice], 578 quick_imshow(stack[z_range[0]:z_range[1],x_range[0]:x_range[1],y_slice],
546 title=title, path=path, save_fig=self.save_figs, save_only=self.save_only, 579 title=title, path=path, save_fig=self.save_figs, save_only=self.save_only,
547 block=self.block) 580 block=self.block)
548 title = f'{basetitle} {res_title} zslice{z_slice}' 581 title = f'{basetitle} {i+1} {res_title} zslice{z_slice}'
549 quick_imshow(stack[z_slice,x_range[0]:x_range[1],y_range[0]:y_range[1]], 582 quick_imshow(stack[z_slice,x_range[0]:x_range[1],y_range[0]:y_range[1]],
550 title=title, path=path, save_fig=self.save_figs, save_only=self.save_only, 583 title=title, path=path, save_fig=self.save_figs, save_only=self.save_only,
551 block=self.block) 584 block=self.block)
552 585
553 # Save test data to file 586 # Save test data to file
738 assert(len(tdf_stack) == 1) # TODO 771 assert(len(tdf_stack) == 1) # TODO
739 tdf_stack = tdf_stack[0] 772 tdf_stack = tdf_stack[0]
740 773
741 # Take median 774 # Take median
742 if tdf_stack.ndim == 2: 775 if tdf_stack.ndim == 2:
743 tdf = tdf_stack 776 tdf = tdf_stack.astype('float64')
744 elif tdf_stack.ndim == 3: 777 elif tdf_stack.ndim == 3:
745 tdf = np.median(tdf_stack, axis=0) 778 tdf = np.median(tdf_stack, axis=0)
746 del tdf_stack 779 del tdf_stack
747 else: 780 else:
748 raise ValueError(f'Invalid tdf_stack shape ({tdf_stack.shape})') 781 raise ValueError(f'Invalid tdf_stack shape ({tdf_stack.shape})')
806 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
807 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.
808 We don’t typically account for them but potentially could. 841 We don’t typically account for them but potentially could.
809 """ 842 """
810 if tbf_stack.ndim == 2: 843 if tbf_stack.ndim == 2:
811 tbf = tbf_stack 844 tbf = tbf_stack.astype('float64')
812 elif tbf_stack.ndim == 3: 845 elif tbf_stack.ndim == 3:
813 tbf = np.median(tbf_stack, axis=0) 846 tbf = np.median(tbf_stack, axis=0)
814 del tbf_stack 847 del tbf_stack
815 else: 848 else:
816 raise ValueError(f'Invalid tbf_stack shape ({tbf_stacks.shape})') 849 raise ValueError(f'Invalid tbf_stack shape ({tbf_stacks.shape})')
821 else: 854 else:
822 logger.warning('Dark field unavailable') 855 logger.warning('Dark field unavailable')
823 856
824 # Set any non-positive values to one 857 # Set any non-positive values to one
825 # (avoid negative bright field values for spikes in dark field) 858 # (avoid negative bright field values for spikes in dark field)
826 tbf[tbf < 1] = 1 859 tbf[tbf < 1.0] = 1.0
827 860
828 # Plot bright field 861 # Plot bright field
829 if self.galaxy_flag: 862 if self.galaxy_flag:
830 quick_imshow(tbf, title='bright field', path='tomo_reduce_plots', 863 quick_imshow(tbf, title='bright field', path='tomo_reduce_plots',
831 save_fig=self.save_figs, save_only=self.save_only) 864 save_fig=self.save_figs, save_only=self.save_only)
871 num_tomo_stacks = len(tomo_fields.scan_numbers) 904 num_tomo_stacks = len(tomo_fields.scan_numbers)
872 theta = tomo_fields.theta_range['start'] 905 theta = tomo_fields.theta_range['start']
873 906
874 # Select image bounds 907 # Select image bounds
875 title = f'tomography image at theta={round(theta, 2)+0}' 908 title = f'tomography image at theta={round(theta, 2)+0}'
876 if (img_x_bounds is not None and not is_index_range(img_x_bounds, ge=0,
877 le=first_image.shape[0])):
878 raise ValueError(f'Invalid parameter img_x_bounds ({img_x_bounds})')
879 if nxentry.instrument.source.attrs['station'] in ('id1a3', 'id3a'): 909 if nxentry.instrument.source.attrs['station'] in ('id1a3', 'id3a'):
880 pixel_size = nxentry.instrument.detector.x_pixel_size 910 pixel_size = nxentry.instrument.detector.x_pixel_size
881 # Try to get a fit from the bright field 911 # Try to get a fit from the bright field
882 tbf = np.asarray(reduced_data.data.bright_field) 912 tbf = np.asarray(reduced_data.data.bright_field)
883 tbf_shape = tbf.shape 913 tbf_shape = tbf.shape
975 x_sum_min = x_sum.min() 1005 x_sum_min = x_sum.min()
976 x_sum_max = x_sum.max() 1006 x_sum_max = x_sum.max()
977 if self.galaxy_flag: 1007 if self.galaxy_flag:
978 if img_x_bounds is None: 1008 if img_x_bounds is None:
979 img_x_bounds = (0, first_image.shape[0]) 1009 img_x_bounds = (0, first_image.shape[0])
1010 elif is_int_pair(img_x_bounds, ge=-1, le=first_image.shape[0]):
1011 img_x_bounds = tuple(img_x_bounds)
1012 if img_x_bounds[0] == -1:
1013 img_x_bounds = (0, img_x_bounds[1])
1014 if img_x_bounds[1] == -1:
1015 img_x_bounds = (img_x_bounds[0], first_image.shape[0])
1016 if not is_index_range(img_x_bounds, ge=0, le=first_image.shape[0]):
1017 raise ValueError(f'Invalid parameter img_x_bounds ({img_x_bounds})')
980 else: 1018 else:
981 quick_imshow(first_image, title=title) 1019 quick_imshow(first_image, title=title)
982 print('Select vertical data reduction range from first tomography image') 1020 print('Select vertical data reduction range from first tomography image')
983 img_x_bounds = select_image_bounds(first_image, 0, title=title) 1021 img_x_bounds = select_image_bounds(first_image, 0, title=title)
984 clear_imshow(title) 1022 clear_imshow(title)
1029 def _gen_tomo(self, nxentry, reduced_data): 1067 def _gen_tomo(self, nxentry, reduced_data):
1030 """Generate tomography fields. 1068 """Generate tomography fields.
1031 """ 1069 """
1032 # Get full bright field 1070 # Get full bright field
1033 tbf = np.asarray(reduced_data.data.bright_field) 1071 tbf = np.asarray(reduced_data.data.bright_field)
1034 tbf_shape = tbf.shape 1072 img_shape = tbf.shape
1035 1073
1036 # Get image bounds 1074 # Get image bounds
1037 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])))
1038 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
1039 1081
1040 # Get resized dark field 1082 # Get resized dark field
1041 # if 'dark_field' in data: 1083 if 'dark_field' in reduced_data.data:
1042 # tbf = np.asarray(reduced_data.data.dark_field[ 1084 if resize_flag:
1043 # 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[
1044 # else: 1086 img_x_bounds[0]:img_x_bounds[1],img_y_bounds[0]:img_y_bounds[1]])
1045 # logger.warning('Dark field unavailable') 1087 else:
1046 # tdf = None 1088 tdf = np.asarray(reduced_data.data.dark_field)
1047 tdf = None 1089 else:
1090 logger.warning('Dark field unavailable')
1091 tdf = None
1048 1092
1049 # Resize bright field 1093 # Resize bright field
1050 if img_x_bounds != (0, tbf.shape[0]) or img_y_bounds != (0, tbf.shape[1]): 1094 if resize_flag:
1051 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]]
1052 1096
1053 # Get the tomography images 1097 # Get the tomography images
1054 image_key = nxentry.instrument.detector.get('image_key', None) 1098 image_key = nxentry.instrument.detector.get('image_key', None)
1055 if image_key and 'data' in nxentry.instrument.detector: 1099 if image_key and 'data' in nxentry.instrument.detector:
1105 if self.galaxy_flag: 1149 if self.galaxy_flag:
1106 path = 'tomo_reduce_plots' 1150 path = 'tomo_reduce_plots'
1107 else: 1151 else:
1108 path = self.output_folder 1152 path = self.output_folder
1109 for i, tomo_stack in enumerate(tomo_stacks): 1153 for i, tomo_stack in enumerate(tomo_stacks):
1110 # Resize the tomography images 1154 # Resize the tomography images as needed
1111 # 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
1112 if img_x_bounds != (0, tbf.shape[0]) or img_y_bounds != (0, tbf.shape[1]): 1156 if resize_flag:
1113 t0 = time() 1157 t0 = time()
1114 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],
1115 img_y_bounds[0]:img_y_bounds[1]].astype('float64') 1159 img_y_bounds[0]:img_y_bounds[1]].astype('float64')
1116 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')
1117 1163
1118 # Subtract dark field 1164 # Subtract dark field
1119 if tdf is not None: 1165 if tdf is not None:
1120 t0 = time() 1166 t0 = time()
1121 with set_numexpr_threads(self.num_core): 1167 with set_numexpr_threads(self.num_core):