Mercurial > repos > rv43 > tomo

comparison tomo.py @ 38:c09b8ee8f68f 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 62ea2c5a3629c5e4b5b59fadc282214326f093aa"
author rv43
date Tue, 19 Apr 2022 18:28:07 +0000
parents 51dc93d0806a
children 8a3036b0a34c
comparison
equal deleted inserted replaced
37:51dc93d0806a 38:c09b8ee8f68f
1145 t0 = time() 1145 t0 = time()
1146 # recon_sinogram = filters.gaussian(recon_sinogram, 3.0) 1146 # recon_sinogram = filters.gaussian(recon_sinogram, 3.0)
1147 recon_sinogram = spi.gaussian_filter(recon_sinogram, 0.5) 1147 recon_sinogram = spi.gaussian_filter(recon_sinogram, 0.5)
1148 recon_clean = np.expand_dims(recon_sinogram, axis=0) 1148 recon_clean = np.expand_dims(recon_sinogram, axis=0)
1149 del recon_sinogram 1149 del recon_sinogram
1150 logging.info(f'tomopy.misc.corr.remove_ring start on {num_core} cores') 1150 t1 = time()
1151 logging.info(f'running tomopy.misc.corr.remove_ring on {num_core} cores ...')
1151 recon_clean = tomopy.misc.corr.remove_ring(recon_clean, rwidth=17, ncore=num_core) 1152 recon_clean = tomopy.misc.corr.remove_ring(recon_clean, rwidth=17, ncore=num_core)
1152 logging.info('tomopy.misc.corr.remove_ring end') 1153 logging.info('... tomopy.misc.corr.remove_ring took {time()-t1:.2f} seconds!')
1153 logging.debug(f'filtering and removing ring artifact took {time()-t0:.2f} seconds!') 1154 logging.debug(f'filtering and removing ring artifact took {time()-t0:.2f} seconds!')
1154 return recon_clean 1155 return recon_clean
1155 1156
1156 def _plotEdgesOnePlane(self, recon_plane, title, path=None, weight=0.001): 1157 def _plotEdgesOnePlane(self, recon_plane, title, path=None, weight=0.001):
1157 # RV parameters for the denoise, gaussian, and ring removal will be different for different feature sizes 1158 # RV parameters for the denoise, gaussian, and ring removal will be different for different feature sizes
1183 # need index order column,theta for iradon, so take transpose 1184 # need index order column,theta for iradon, so take transpose
1184 sinogram_T = sinogram.T 1185 sinogram_T = sinogram.T
1185 center = sinogram.shape[1]/2 1186 center = sinogram.shape[1]/2
1186 1187
1187 # try automatic center finding routines for initial value 1188 # try automatic center finding routines for initial value
1189 t0 = time()
1188 if num_core > 24: 1190 if num_core > 24:
1189 logging.info('tomopy.find_center_vo start on 24 cores') 1191 logging.info('running tomopy.find_center_vo on 24 cores ...')
1190 tomo_center = tomopy.find_center_vo(sinogram, ncore=24) 1192 tomo_center = tomopy.find_center_vo(sinogram, ncore=24)
1191 else: 1193 else:
1192 logging.info(f'tomopy.find_center_vo start on {num_core} cores') 1194 logging.info(f'running tomopy.find_center_vo on {num_core} cores ...')
1193 tomo_center = tomopy.find_center_vo(sinogram, ncore=num_core) 1195 tomo_center = tomopy.find_center_vo(sinogram, ncore=num_core)
1194 logging.info('tomopy.find_center_vo end') 1196 logging.info('... tomopy.find_center_vo took {time()-t0:.2f} seconds!')
1195 center_offset_vo = tomo_center-center 1197 center_offset_vo = tomo_center-center
1196 if self.test_mode: 1198 if self.test_mode:
1197 logging.info(f'Center at row {row} using Nghia Vo’s method = {center_offset_vo:.2f}') 1199 logging.info(f'Center at row {row} using Nghia Vo’s method = {center_offset_vo:.2f}')
1198 del sinogram_T 1200 del sinogram_T
1199 return float(center_offset_vo) 1201 return float(center_offset_vo)
1200 elif self.galaxy_flag: 1202 elif self.galaxy_flag:
1201 logging.info(f'Center at row {row} using Nghia Vo’s method = {center_offset_vo:.2f}') 1203 logging.info(f'Center at row {row} using Nghia Vo’s method = {center_offset_vo:.2f}')
1204 t0 = time()
1205 logging.info('running self._reconstructOnePlane on {num_core} cores ...')
1202 recon_plane = self._reconstructOnePlane(sinogram_T, tomo_center, thetas_deg, 1206 recon_plane = self._reconstructOnePlane(sinogram_T, tomo_center, thetas_deg,
1203 eff_pixel_size, cross_sectional_dim, False, num_core) 1207 eff_pixel_size, cross_sectional_dim, False, num_core)
1208 logging.info('... self._reconstructOnePlane took {time()-t0:.2f} seconds!')
1204 title = f'edges row{row} center offset{center_offset_vo:.2f} Vo' 1209 title = f'edges row{row} center offset{center_offset_vo:.2f} Vo'
1205 self._plotEdgesOnePlane(recon_plane, title, path='find_center_pngs') 1210 self._plotEdgesOnePlane(recon_plane, title, path='find_center_pngs')
1206 del recon_plane 1211 del recon_plane
1207 if not galaxy_param['center_type_selector']: 1212 if not galaxy_param['center_type_selector']:
1208 del sinogram_T 1213 del sinogram_T
1273 center_offsets = np.linspace(center_offset_low, center_offset_upp, num_center_offset) 1278 center_offsets = np.linspace(center_offset_low, center_offset_upp, num_center_offset)
1274 for center_offset in center_offsets: 1279 for center_offset in center_offsets:
1275 if center_offset == center_offset_vo: 1280 if center_offset == center_offset_vo:
1276 continue 1281 continue
1277 logging.info(f'center_offset = {center_offset:.2f}') 1282 logging.info(f'center_offset = {center_offset:.2f}')
1283 t0 = time()
1284 logging.info('running self._reconstructOnePlane on {num_core} cores ...')
1278 recon_plane = self._reconstructOnePlane(sinogram_T, center_offset+center, 1285 recon_plane = self._reconstructOnePlane(sinogram_T, center_offset+center,
1279 thetas_deg, eff_pixel_size, cross_sectional_dim, False, num_core) 1286 thetas_deg, eff_pixel_size, cross_sectional_dim, False, num_core)
1287 logging.info('... self._reconstructOnePlane took {time()-t0:.2f} seconds!')
1280 title = f'edges row{row} center_offset{center_offset:.2f}' 1288 title = f'edges row{row} center_offset{center_offset:.2f}'
1281 if self.galaxy_flag: 1289 if self.galaxy_flag:
1282 self._plotEdgesOnePlane(recon_plane, title, path='find_center_pngs') 1290 self._plotEdgesOnePlane(recon_plane, title, path='find_center_pngs')
1283 else: 1291 else:
1284 self._plotEdgesOnePlane(recon_plane, title) 1292 self._plotEdgesOnePlane(recon_plane, title)