comparison tomo.py @ 15:1bcca1f2adb4 draft

"planemo upload for repository https://github.com/rolfverberg/galaxytools commit 38c19bf5addbc46f45d598f981bb1a48f7bca691"
author rv43
date Wed, 13 Apr 2022 16:02:44 +0000
parents 40395e60d2be
children 7f723407beb3
comparison
equal deleted inserted replaced
14:50c8d19d5f89 15:1bcca1f2adb4
981 t0 = time() 981 t0 = time()
982 logging.info(f'Saving {tomo_file} ...') 982 logging.info(f'Saving {tomo_file} ...')
983 np.save(tomo_file, stack) 983 np.save(tomo_file, stack)
984 logging.info(f'... done in {time()-t0:.2f} seconds!') 984 logging.info(f'... done in {time()-t0:.2f} seconds!')
985 985
986 def _genTomo(self, tomo_stack_files, available_stacks): 986 def _genTomo(self, tomo_stack_files, available_stacks, num_core=None):
987 """Generate tomography fields. 987 """Generate tomography fields.
988 """ 988 """
989 if num_core is None:
990 num_core = self.num_core
989 stacks = self.config['stack_info']['stacks'] 991 stacks = self.config['stack_info']['stacks']
990 assert(len(self.tomo_stacks) == self.config['stack_info']['num']) 992 assert(len(self.tomo_stacks) == self.config['stack_info']['num'])
991 assert(len(self.tomo_stacks) == len(stacks)) 993 assert(len(self.tomo_stacks) == len(stacks))
992 if len(available_stacks) != len(stacks): 994 if len(available_stacks) != len(stacks):
993 logging.warning('Illegal dimension of available_stacks in _genTomo'+ 995 logging.warning('Illegal dimension of available_stacks in _genTomo'+
1443 self.tomo_stacks[i] = f[f'set_{stack["index"]}'] 1445 self.tomo_stacks[i] = f[f'set_{stack["index"]}']
1444 logging.info(f'loaded stack {i}: index = {stack["index"]}, shape = '+ 1446 logging.info(f'loaded stack {i}: index = {stack["index"]}, shape = '+
1445 f'{self.tomo_stacks[i].shape}') 1447 f'{self.tomo_stacks[i].shape}')
1446 logging.info(f'... done in {time()-t0:.2f} seconds!') 1448 logging.info(f'... done in {time()-t0:.2f} seconds!')
1447 1449
1448 def genTomoStacks(self, tdf_files=None, tbf_files=None, tomo_stack_files=None, 1450 def genTomoStacks(self, galaxy_param=None, num_core=None):
1449 dark_field_pngname=None, bright_field_pngname=None, tomo_field_pngname=None,
1450 detectorbounds_pngname=None, output_name=None):
1451 """Preprocess tomography images. 1451 """Preprocess tomography images.
1452 """ 1452 """
1453 if num_core is None:
1454 num_core = self.num_core
1453 # Try loading any already preprocessed stacks (skip in Galaxy) 1455 # Try loading any already preprocessed stacks (skip in Galaxy)
1454 # preprocessed stack order for each one in stack: row,theta,column 1456 # preprocessed stack order for each one in stack: row,theta,column
1455 stack_info = self.config['stack_info'] 1457 stack_info = self.config['stack_info']
1456 stacks = stack_info['stacks'] 1458 stacks = stack_info['stacks']
1457 num_tomo_stacks = stack_info['num'] 1459 num_tomo_stacks = stack_info['num']
1458 assert(num_tomo_stacks == len(self.tomo_stacks)) 1460 assert(num_tomo_stacks == len(self.tomo_stacks))
1459 available_stacks = [False]*num_tomo_stacks 1461 available_stacks = [False]*num_tomo_stacks
1460 if self.galaxy_flag: 1462 if self.galaxy_flag:
1461 assert(tdf_files is None or isinstance(tdf_files, list)) 1463 assert(isinstance(galaxy_param, dict))
1462 assert(isinstance(tbf_files, list)) 1464 tdf_files = galaxy_param['tdf_files']
1463 assert(isinstance(tomo_stack_files, list)) 1465 tbf_files = galaxy_param['tbf_files']
1466 tomo_stack_files = galaxy_param['tomo_stack_files']
1464 assert(num_tomo_stacks == len(tomo_stack_files)) 1467 assert(num_tomo_stacks == len(tomo_stack_files))
1465 assert(isinstance(dark_field_pngname, str)) 1468 else:
1466 assert(isinstance(bright_field_pngname, str)) 1469 if galaxy_param:
1467 assert(isinstance(tomo_field_pngname, str)) 1470 logging.warning('Ignoring galaxy_param in findCenters (only for Galaxy)')
1468 assert(isinstance(detectorbounds_pngname, str)) 1471 galaxy_param = None
1469 assert(isinstance(output_name, str))
1470 else:
1471 if tdf_files:
1472 logging.warning('Ignoring tdf_files in genTomoStacks (only for Galaxy)')
1473 if tbf_files:
1474 logging.warning('Ignoring tbf_files in genTomoStacks (only for Galaxy)')
1475 if tomo_stack_files:
1476 logging.warning('Ignoring tomo_stack_files in genTomoStacks (only for Galaxy)')
1477 tdf_files, tbf_files, tomo_stack_files = self.findImageFiles() 1472 tdf_files, tbf_files, tomo_stack_files = self.findImageFiles()
1478 if not self.is_valid: 1473 if not self.is_valid:
1479 return 1474 return
1480 for i,stack in enumerate(stacks): 1475 for i,stack in enumerate(stacks):
1481 if not self.tomo_stacks[i].size and stack.get('preprocessed', False): 1476 if not self.tomo_stacks[i].size and stack.get('preprocessed', False):
1482 self.tomo_stacks[i], available_stacks[i] = \ 1477 self.tomo_stacks[i], available_stacks[i] = \
1483 self._loadTomo('red stack', stack['index']) 1478 self._loadTomo('red stack', stack['index'])
1484 if dark_field_pngname:
1485 logging.warning('Ignoring dark_field_pngname in genTomoStacks (only for Galaxy)')
1486 if bright_field_pngname:
1487 logging.warning('Ignoring bright_field_pngname in genTomoStacks (only for Galaxy)')
1488 if tomo_field_pngname:
1489 logging.warning('Ignoring tomo_field_pngname in genTomoStacks (only for Galaxy)')
1490 if detectorbounds_pngname:
1491 logging.warning('Ignoring detectorbounds_pngname in genTomoStacks '+
1492 '(only used in Galaxy)')
1493 if output_name:
1494 logging.warning('Ignoring output_name in genTomoStacks '+
1495 '(only used in Galaxy)')
1496 1479
1497 # Preprocess any unloaded stacks 1480 # Preprocess any unloaded stacks
1498 if False in available_stacks: 1481 if False in available_stacks:
1499 logging.debug('Preprocessing tomography images') 1482 logging.debug('Preprocessing tomography images')
1500 1483
1504 if not self.is_valid: 1487 if not self.is_valid:
1505 return 1488 return
1506 1489
1507 # Generate dark field 1490 # Generate dark field
1508 if tdf_files: 1491 if tdf_files:
1509 self._genDark(tdf_files, dark_field_pngname) 1492 self._genDark(tdf_files, galaxy_param['dark_field_pngname'])
1510 1493
1511 # Generate bright field 1494 # Generate bright field
1512 self._genBright(tbf_files, bright_field_pngname) 1495 self._genBright(tbf_files, galaxy_param['bright_field_pngname'])
1513 1496
1514 # Set vertical detector bounds for image stack 1497 # Set vertical detector bounds for image stack
1515 self._setDetectorBounds(tomo_stack_files, tomo_field_pngname, detectorbounds_pngname) 1498 self._setDetectorBounds(tomo_stack_files, galaxy_param['tomo_field_pngname'],
1499 galaxy_param['detectorbounds_pngname'])
1516 1500
1517 # Set zoom and/or theta skip to reduce memory the requirement 1501 # Set zoom and/or theta skip to reduce memory the requirement
1518 self._setZoomOrSkip() 1502 self._setZoomOrSkip()
1519 1503
1520 # Generate tomography fields 1504 # Generate tomography fields
1521 self._genTomo(tomo_stack_files, available_stacks) 1505 self._genTomo(tomo_stack_files, available_stacks, num_core)
1522 1506
1523 # Save tomography stack to file 1507 # Save tomography stack to file
1524 if self.galaxy_flag: 1508 if self.galaxy_flag:
1525 t0 = time() 1509 t0 = time()
1510 output_name = galaxy_param['output_name']
1526 logging.info(f'Saving preprocessed tomography stack to {output_name} ...') 1511 logging.info(f'Saving preprocessed tomography stack to {output_name} ...')
1527 save_stacks = {f'set_{stack["index"]}':tomo_stack 1512 save_stacks = {f'set_{stack["index"]}':tomo_stack
1528 for stack,tomo_stack in zip(stacks,self.tomo_stacks)} 1513 for stack,tomo_stack in zip(stacks,self.tomo_stacks)}
1529 np.savez(output_name, **save_stacks) 1514 np.savez(output_name, **save_stacks)
1530 logging.info(f'... done in {time()-t0:.2f} seconds!') 1515 logging.info(f'... done in {time()-t0:.2f} seconds!')
1887 del plane1, plane2, plane1_shifted 1872 del plane1, plane2, plane1_shifted
1888 1873
1889 # Update config file 1874 # Update config file
1890 self.config = msnc.update('config.txt', 'check_centers', True, 'find_centers') 1875 self.config = msnc.update('config.txt', 'check_centers', True, 'find_centers')
1891 1876
1892 def reconstructTomoStacks(self, output_name=None, num_core=None): 1877 def reconstructTomoStacks(self, galaxy_param=None, num_core=None):
1893 """Reconstruct tomography stacks. 1878 """Reconstruct tomography stacks.
1894 """ 1879 """
1895 if num_core is None: 1880 if num_core is None:
1896 num_core = self.num_core 1881 num_core = self.num_core
1897 logging.debug('Reconstruct tomography stacks') 1882 logging.debug('Reconstruct tomography stacks')
1898 stacks = self.config['stack_info']['stacks'] 1883 stacks = self.config['stack_info']['stacks']
1899 assert(len(self.tomo_stacks) == self.config['stack_info']['num']) 1884 assert(len(self.tomo_stacks) == self.config['stack_info']['num'])
1900 assert(len(self.tomo_stacks) == len(stacks)) 1885 assert(len(self.tomo_stacks) == len(stacks))
1901 assert(len(self.tomo_recon_stacks) == len(stacks)) 1886 assert(len(self.tomo_recon_stacks) == len(stacks))
1902 if self.galaxy_flag: 1887 if self.galaxy_flag:
1903 assert(isinstance(output_name, str)) 1888 assert(isinstance(galaxy_param, dict))
1904 else: 1889 # Get rotation axis centers
1905 if output_name: 1890 center_offsets = galaxy_param['center_offsets']
1906 logging.warning('Ignoring output_name in reconstructTomoStacks '+ 1891 assert(isinstance(center_offsets, list) and len(center_offsets) == 2)
1907 '(only used in Galaxy)') 1892 lower_center_offset = center_offsets[0]
1893 assert(msnc.is_num(lower_center_offset))
1894 upper_center_offset = center_offsets[1]
1895 assert(msnc.is_num(upper_center_offset))
1896 else:
1897 if galaxy_param:
1898 logging.warning('Ignoring galaxy_param in reconstructTomoStacks (only for Galaxy)')
1899 galaxy_param = None
1900 lower_center_offset = None
1901 upper_center_offset = None
1908 1902
1909 # Get rotation axis rows and centers 1903 # Get rotation axis rows and centers
1910 find_center = self.config['find_center'] 1904 find_center = self.config['find_center']
1911 lower_row = find_center.get('lower_row') 1905 lower_row = find_center.get('lower_row')
1912 if lower_row is None: 1906 if lower_row is None:
1913 logging.error('Unable to read lower_row from config') 1907 logging.error('Unable to read lower_row from config')
1914 return 1908 return
1915 lower_center_offset = find_center.get('lower_center_offset')
1916 if lower_center_offset is None:
1917 logging.error('Unable to read lower_center_offset from config')
1918 return
1919 upper_row = find_center.get('upper_row') 1909 upper_row = find_center.get('upper_row')
1920 if upper_row is None: 1910 if upper_row is None:
1921 logging.error('Unable to read upper_row from config') 1911 logging.error('Unable to read upper_row from config')
1922 return 1912 return
1923 upper_center_offset = find_center.get('upper_center_offset')
1924 if upper_center_offset is None:
1925 logging.error('Unable to read upper_center_offset from config')
1926 return
1927 logging.debug(f'lower_row = {lower_row} upper_row = {upper_row}') 1913 logging.debug(f'lower_row = {lower_row} upper_row = {upper_row}')
1928 assert(lower_row < upper_row) 1914 assert(lower_row < upper_row)
1915 if lower_center_offset is None:
1916 lower_center_offset = find_center.get('lower_center_offset')
1917 if lower_center_offset is None:
1918 logging.error('Unable to read lower_center_offset from config')
1919 return
1920 if upper_center_offset is None:
1921 upper_center_offset = find_center.get('upper_center_offset')
1922 if upper_center_offset is None:
1923 logging.error('Unable to read upper_center_offset from config')
1924 return
1929 center_slope = (upper_center_offset-lower_center_offset)/(upper_row-lower_row) 1925 center_slope = (upper_center_offset-lower_center_offset)/(upper_row-lower_row)
1930 1926
1931 # Set thetas (in radians) 1927 # Set thetas (in radians)
1932 theta_range = self.config['theta_range'] 1928 theta_range = self.config['theta_range']
1933 theta_start = theta_range['start'] 1929 theta_start = theta_range['start']
2002 self.cf.saveFile(self.config_out) 1998 self.cf.saveFile(self.config_out)
2003 1999
2004 # Save reconstructed tomography stack to file 2000 # Save reconstructed tomography stack to file
2005 if self.galaxy_flag: 2001 if self.galaxy_flag:
2006 t0 = time() 2002 t0 = time()
2003 output_name = galaxy_param['output_name']
2007 logging.info(f'Saving reconstructed tomography stack to {output_name} ...') 2004 logging.info(f'Saving reconstructed tomography stack to {output_name} ...')
2008 save_stacks = {f'set_{stack["index"]}':tomo_stack 2005 save_stacks = {f'set_{stack["index"]}':tomo_stack
2009 for stack,tomo_stack in zip(stacks,self.tomo_recon_stacks)} 2006 for stack,tomo_stack in zip(stacks,self.tomo_recon_stacks)}
2010 np.savez(output_name, **save_stacks) 2007 np.savez(output_name, **save_stacks)
2011 logging.info(f'... done in {time()-t0:.2f} seconds!') 2008 logging.info(f'... done in {time()-t0:.2f} seconds!')