comparison tomo.py @ 40:fa94fe25ca46 draft

"planemo upload for repository https://github.com/rolfverberg/galaxytools commit 9978af8877d8e9722d252ab9bee4cbd4654aa100"
author rv43
date Tue, 19 Apr 2022 19:11:59 +0000
parents 8a3036b0a34c
children ef5c2f7b49ec
comparison
equal deleted inserted replaced
39:8a3036b0a34c 40:fa94fe25ca46
913 else: 913 else:
914 num_theta_skip = preprocess['num_theta_skip'] 914 num_theta_skip = preprocess['num_theta_skip']
915 if not msnc.is_int(num_theta_skip, 0): 915 if not msnc.is_int(num_theta_skip, 0):
916 msnc.illegal_value('num_theta_skip', num_theta_skip, 'config file') 916 msnc.illegal_value('num_theta_skip', num_theta_skip, 'config file')
917 num_theta_skip = 0 917 num_theta_skip = 0
918 logging.info(f'zoom_perc = {zoom_perc}') 918 logging.debug(f'zoom_perc = {zoom_perc}')
919 logging.info(f'num_theta_skip = {num_theta_skip}') 919 logging.debug(f'num_theta_skip = {num_theta_skip}')
920 920
921 # Update config and save to file 921 # Update config and save to file
922 if preprocess is None: 922 if preprocess is None:
923 self.cf.config['preprocess'] = {'zoom_perc' : zoom_perc, 923 self.cf.config['preprocess'] = {'zoom_perc' : zoom_perc,
924 'num_theta_skip' : num_theta_skip} 924 'num_theta_skip' : num_theta_skip}
1331 else: 1331 else:
1332 if center_offsets.size != row_bounds[1]-row_bounds[0]: 1332 if center_offsets.size != row_bounds[1]-row_bounds[0]:
1333 raise ValueError('center_offsets dimension mismatch in reconstructOneTomoStack') 1333 raise ValueError('center_offsets dimension mismatch in reconstructOneTomoStack')
1334 centers = center_offsets 1334 centers = center_offsets
1335 centers += tomo_stack.shape[2]/2 1335 centers += tomo_stack.shape[2]/2
1336 # RV hangs here with more than 24 cores and sge_64G_4
1337 if True: 1336 if True:
1338 t0 = time() 1337 t0 = time()
1339 logging.info(f'running tomopy.prep.stripe.remove_stripe_fw on {num_core} cores ...') 1338 if num_core > 24:
1340 tomo_stack = tomopy.prep.stripe.remove_stripe_fw( 1339 logging.info(f'running tomopy.prep.stripe.remove_stripe_fw on 24 cores ...')
1341 tomo_stack[row_bounds[0]:row_bounds[1]], sigma=sigma, ncore=num_core) 1340 tomo_stack = tomopy.prep.stripe.remove_stripe_fw(
1341 tomo_stack[row_bounds[0]:row_bounds[1]], sigma=sigma, ncore=24)
1342 else:
1343 logging.info(f'running tomopy.prep.stripe.remove_stripe_fw on {num_core} cores ...')
1344 tomo_stack = tomopy.prep.stripe.remove_stripe_fw(
1345 tomo_stack[row_bounds[0]:row_bounds[1]], sigma=sigma, ncore=num_core)
1342 logging.info(f'... tomopy.prep.stripe.remove_stripe_fw took {time()-t0:.2f} seconds!') 1346 logging.info(f'... tomopy.prep.stripe.remove_stripe_fw took {time()-t0:.2f} seconds!')
1343 else: 1347 else:
1344 tomo_stack = tomo_stack[row_bounds[0]:row_bounds[1]] 1348 tomo_stack = tomo_stack[row_bounds[0]:row_bounds[1]]
1345 logging.info('performing initial reconstruction') 1349 logging.info('performing initial reconstruction')
1346 t0 = time() 1350 t0 = time()
1565 def findCenters(self, galaxy_param=None, num_core=None): 1569 def findCenters(self, galaxy_param=None, num_core=None):
1566 """Find rotation axis centers for the tomography stacks. 1570 """Find rotation axis centers for the tomography stacks.
1567 """ 1571 """
1568 if num_core is None: 1572 if num_core is None:
1569 num_core = self.num_core 1573 num_core = self.num_core
1570 # logging.info(f'num_core available = {num_core}') 1574 logging.info(f'num_core = {num_core}')
1571 # if num_core > 24:
1572 # num_core = 24
1573 logging.info(f'num_core used = {num_core}')
1574 logging.debug('Find centers for tomography stacks') 1575 logging.debug('Find centers for tomography stacks')
1575 stacks = self.config['stack_info']['stacks'] 1576 stacks = self.config['stack_info']['stacks']
1576 available_stacks = [stack['index'] for stack in stacks if stack.get('preprocessed', False)] 1577 available_stacks = [stack['index'] for stack in stacks if stack.get('preprocessed', False)]
1577 logging.debug('Available stacks: {available_stacks}') 1578 logging.debug('Available stacks: {available_stacks}')
1578 if self.galaxy_flag: 1579 if self.galaxy_flag:
1590 set_center = galaxy_param['set_center'] 1591 set_center = galaxy_param['set_center']
1591 else: 1592 else:
1592 logging.error('Illegal center_type_selector entry in galaxy_param '+ 1593 logging.error('Illegal center_type_selector entry in galaxy_param '+
1593 f'({center_type_selector})') 1594 f'({center_type_selector})')
1594 galaxy_param['center_type_selector'] = None 1595 galaxy_param['center_type_selector'] = None
1595 logging.info(f'row_bounds = {row_bounds}') 1596 logging.debug(f'row_bounds = {row_bounds}')
1596 logging.info(f'center_rows = {center_rows}') 1597 logging.debug(f'center_rows = {center_rows}')
1597 logging.info(f'center_type_selector = {center_type_selector}') 1598 logging.debug(f'center_type_selector = {center_type_selector}')
1598 else: 1599 else:
1599 if galaxy_param: 1600 if galaxy_param:
1600 logging.warning('Ignoring galaxy_param in findCenters (only for Galaxy)') 1601 logging.warning('Ignoring galaxy_param in findCenters (only for Galaxy)')
1601 galaxy_param = None 1602 galaxy_param = None
1602 1603
1720 logging.error('CHECK/FIX FOR GALAXY') 1721 logging.error('CHECK/FIX FOR GALAXY')
1721 tomo_ref_heights = [stack['ref_height'] for stack in stacks] 1722 tomo_ref_heights = [stack['ref_height'] for stack in stacks]
1722 n1 = int((1.+(tomo_ref_heights[0]+center_stack.shape[0]*eff_pixel_size- 1723 n1 = int((1.+(tomo_ref_heights[0]+center_stack.shape[0]*eff_pixel_size-
1723 tomo_ref_heights[1])/eff_pixel_size)/2) 1724 tomo_ref_heights[1])/eff_pixel_size)/2)
1724 n2 = center_stack.shape[0]-n1 1725 n2 = center_stack.shape[0]-n1
1725 logging.info(f'n1 = {n1}, n2 = {n2} (n2-n1) = {(n2-n1)*eff_pixel_size:.3f} mm') 1726 logging.debug(f'n1 = {n1}, n2 = {n2} (n2-n1) = {(n2-n1)*eff_pixel_size:.3f} mm')
1726 if not center_stack.size: 1727 if not center_stack.size:
1727 RuntimeError('Center stack not loaded') 1728 RuntimeError('Center stack not loaded')
1728 if not self.test_mode and not self.galaxy_flag: 1729 if not self.test_mode and not self.galaxy_flag:
1729 tmp = center_stack[:,0,:] 1730 tmp = center_stack[:,0,:]
1730 tmp_max = tmp.max() 1731 tmp_max = tmp.max()
1768 if use_row != 'no': 1769 if use_row != 'no':
1769 center_offset = find_center.get('lower_center_offset') 1770 center_offset = find_center.get('lower_center_offset')
1770 if msnc.is_num(center_offset): 1771 if msnc.is_num(center_offset):
1771 use_center = pyip.inputYesNo('Current lower center offset = '+ 1772 use_center = pyip.inputYesNo('Current lower center offset = '+
1772 f'{center_offset}, use this value ([y]/n)? ', blank=True) 1773 f'{center_offset}, use this value ([y]/n)? ', blank=True)
1773 #logging.info(f'use_center = {use_center}')
1774 #logging.info(f'use_row = {use_row}')
1775 if use_center == 'no': 1774 if use_center == 'no':
1776 if use_row == 'no': 1775 if use_row == 'no':
1777 if not self.test_mode: 1776 if not self.test_mode:
1778 msnc.quickImshow(center_stack[:,0,:], title=f'theta={theta_start}', 1777 msnc.quickImshow(center_stack[:,0,:], title=f'theta={theta_start}',
1779 aspect='auto') 1778 aspect='auto')
1799 # Upper row center 1798 # Upper row center
1800 use_row = 'no' 1799 use_row = 'no'
1801 use_center = 'no' 1800 use_center = 'no'
1802 row = center_rows[1] 1801 row = center_rows[1]
1803 if self.test_mode or self.galaxy_flag: 1802 if self.test_mode or self.galaxy_flag:
1804 logging.info(f'row = {row} lower_row = {lower_row} n2 = {n2}')
1805 assert(msnc.is_int(row, lower_row+1, n2-1)) 1803 assert(msnc.is_int(row, lower_row+1, n2-1))
1806 if msnc.is_int(row, lower_row+1, n2-1): 1804 if msnc.is_int(row, lower_row+1, n2-1):
1807 if self.test_mode or self.galaxy_flag: 1805 if self.test_mode or self.galaxy_flag:
1808 use_row = 'yes' 1806 use_row = 'yes'
1809 else: 1807 else:
1815 if use_row != 'no': 1813 if use_row != 'no':
1816 center_offset = find_center.get('upper_center_offset') 1814 center_offset = find_center.get('upper_center_offset')
1817 if msnc.is_num(center_offset): 1815 if msnc.is_num(center_offset):
1818 use_center = pyip.inputYesNo('Current upper center offset = '+ 1816 use_center = pyip.inputYesNo('Current upper center offset = '+
1819 f'{center_offset}, use this value ([y]/n)? ', blank=True) 1817 f'{center_offset}, use this value ([y]/n)? ', blank=True)
1820 #logging.info(f'use_center = {use_center}')
1821 #logging.info(f'use_row = {use_row}')
1822 if use_center == 'no': 1818 if use_center == 'no':
1823 if use_row == 'no': 1819 if use_row == 'no':
1824 if not self.test_mode: 1820 if not self.test_mode:
1825 msnc.quickImshow(center_stack[:,0,:], title=f'theta={theta_start}', 1821 msnc.quickImshow(center_stack[:,0,:], title=f'theta={theta_start}',
1826 aspect='auto') 1822 aspect='auto')