From be242888cd554002870fa5b8c3e716a8c86bcff4 Mon Sep 17 00:00:00 2001 From: procyontao Date: Sat, 25 Sep 2021 12:47:17 +0800 Subject: [PATCH] preserve pixelsize during deconvolution and prediction --- bin/isonet.py | 111 +++++++++++++++++++++--------------------- bin/predict.py | 14 +++--- util/deconvolution.py | 43 +++++++++------- 3 files changed, 88 insertions(+), 80 deletions(-) diff --git a/bin/isonet.py b/bin/isonet.py index 1c8aab8..6cbab06 100755 --- a/bin/isonet.py +++ b/bin/isonet.py @@ -21,27 +21,26 @@ class ISONET: """ #log_file = "log.txt" - + def prepare_star(self,folder_name, output_star='tomograms.star',pixel_size = 10.0, defocus = 0.0, number_subtomos = 100): """ \nThis command generates a tomograms.star file from a folder containing only tomogram files (.mrc or .rec).\n isonet.py prepare_star folder_name [--output_star] [--pixel_size] [--defocus] [--number_subtomos] :param folder_name: (None) directory containing tomogram(s). Usually 1-5 tomograms are sufficient. :param output_star: (tomograms.star) star file similar to that from "relion". You can modify this file manually or with gui. - :param pixel_size: (10) pixel size in anstroms. Usually you want to bin your tomograms to about 10A pixel size. + :param pixel_size: (10) pixel size in anstroms. Usually you want to bin your tomograms to about 10A pixel size. Too large or too small pixel sizes are not recommanded, since the target resolution on Z-axis of corrected tomograms should be about 30A. - :param defocus: (0.0) defocus in Angstrom. Only need for ctf deconvolution. For phase plate data, you can leave defocus 0. + :param defocus: (0.0) defocus in Angstrom. Only need for ctf deconvolution. For phase plate data, you can leave defocus 0. If you have multiple tomograms with different defocus, please modify them in star file or with gui. - :param number_subtomos: (100) Number of subtomograms to be extracted in later processes. + :param number_subtomos: (100) Number of subtomograms to be extracted in later processes. If you want to extract different number of subtomograms in different tomograms, you can modify them in the star file generated with this command or with gui. - """ + """ md = MetaData() md.addLabels('rlnIndex','rlnMicrographName','rlnPixelSize','rlnDefocus','rlnNumberSubtomo') tomo_list = sorted(os.listdir(folder_name)) i = 0 for tomo in tomo_list: - print(tomo[-4:]) if tomo[-4:] == '.rec' or tomo[-4:] == '.mrc': i+=1 it = Item() @@ -57,13 +56,13 @@ def prepare_subtomo_star(self, folder_name, output_star='subtomo.star', pixel_si """ \nThis command generates a subtomo star file from a folder containing only subtomogram files (.mrc). This command is usually not necessary in the traditional workflow, because "isonet.py extract" will generate this subtomo.star for you.\n - isonet.py prepare_subtomo_star folder_name [--output_star] [--cube_size] + isonet.py prepare_subtomo_star folder_name [--output_star] [--cube_size] :param folder_name: (None) directory containing subtomogram(s). :param output_star: (subtomo.star) output star file for subtomograms, will be used as input in refinement. - :param pixel_size: (10) The pixel size in angstrom of your subtomograms. - :param cube_size: (None) This is the size of the cubic volumes used for training. This values should be smaller than the size of subtomogram. + :param pixel_size: (10) The pixel size in angstrom of your subtomograms. + :param cube_size: (None) This is the size of the cubic volumes used for training. This values should be smaller than the size of subtomogram. And the cube_size should be divisible by 8. If this value isn't set, cube_size is automatically determined as int(subtomo_size / 1.5 + 1)//16 * 16 - """ + """ #TODO check folder valid, logging if not os.path.isdir(folder_name): print("the folder does not exist") @@ -73,7 +72,7 @@ def prepare_subtomo_star(self, folder_name, output_star='subtomo.star', pixel_si subtomo_list = sorted(os.listdir(folder_name)) for i,subtomo in enumerate(subtomo_list): subtomo_name = os.path.join(folder_name,subtomo) - try: + try: with mrcfile.open(subtomo_name, mode='r') as s: crop_size = s.header.nx except: @@ -97,10 +96,10 @@ def prepare_subtomo_star(self, folder_name, output_star='subtomo.star', pixel_si # f.write(str(i+1)+' ' + os.path.join(folder_name,tomo) + '\n') md.write(output_star) - def deconv(self, star_file: str, - deconv_folder:str="./deconv", - snrfalloff: float=None, - deconvstrength: float=None, + def deconv(self, star_file: str, + deconv_folder:str="./deconv", + snrfalloff: float=None, + deconvstrength: float=None, highpassnyquist: float=0.02, tile: tuple=(1,4,4), overlap_rate = 0.25, @@ -109,28 +108,28 @@ def deconv(self, star_file: str, """ \nCTF deconvolution for the tomograms.\n isonet.py deconv star_file [--deconv_folder] [--snrfalloff] [--deconvstrength] [--highpassnyquist] [--tile] [--overlap_rate] [--ncpu] [--tomo_idx] - This step is recommanded because it enhances low resolution information for a better contrast. No need to do deconvolution for phase plate data. + This step is recommanded because it enhances low resolution information for a better contrast. No need to do deconvolution for phase plate data. :param deconv_folder: (./deconv) Folder created to save deconvoluted tomograms. :param star_file: (None) Star file for tomograms. - :param snrfalloff: (1.0) SNR fall rate with the frequency. High values means losing more high frequency. - If this value is not set, the program will look for the parameter in the star file. + :param snrfalloff: (1.0) SNR fall rate with the frequency. High values means losing more high frequency. + If this value is not set, the program will look for the parameter in the star file. If this value is not set and not found in star file, the default value 1.0 will be used. - :param deconvstrength: (1.0) Strength of the deconvolution. - If this value is not set, the program will look for the parameter in the star file. + :param deconvstrength: (1.0) Strength of the deconvolution. + If this value is not set, the program will look for the parameter in the star file. If this value is not set and not found in star file, the default value 1.0 will be used. :param highpassnyquist: (0.02) Highpass filter for at very low frequency. We suggest to keep this default value. :param tile: (1,4,4) The program crop the tomogram in multiple tiles (z,y,x) for multiprocessing and assembly them into one. e.g. (1,2,2) :param overlap_rate: (None) The overlapping rate for adjecent tiles. - :param ncpu: (4) Number of cpus to use. - :param tomo_idx: (None) If this value is set, process only the tomograms listed in this index. e.g. 1,2,4 or 5-10,15,16 - """ + :param ncpu: (4) Number of cpus to use. + :param tomo_idx: (None) If this value is set, process only the tomograms listed in this index. e.g. 1,2,4 or 5-10,15,16 + """ from IsoNet.util.deconvolution import deconv_one logging.basicConfig(format='%(asctime)s, %(levelname)-8s %(message)s', datefmt="%m-%d %H:%M:%S",level=logging.INFO,handlers=[logging.StreamHandler(sys.stdout)]) logging.info('\n######Isonet starts ctf deconvolve######\n') - try: + try: md = MetaData() md.read(star_file) if not 'rlnSnrFalloff' in md.getLabels(): @@ -139,10 +138,10 @@ def deconv(self, star_file: str, md._setItemValue(it,Label('rlnSnrFalloff'),1.0) md._setItemValue(it,Label('rlnDeconvStrength'),1.0) md._setItemValue(it,Label('rlnDeconvTomoName'),None) - + if not os.path.isdir(deconv_folder): os.mkdir(deconv_folder) - + tomo_idx = idx2list(tomo_idx) for it in md: if tomo_idx is None or str(it.rlnIndex) in tomo_idx: @@ -150,16 +149,16 @@ def deconv(self, star_file: str, md._setItemValue(it,Label('rlnSnrFalloff'), snrfalloff) if deconvstrength is not None: md._setItemValue(it,Label('rlnDeconvStrength'),deconvstrength) - + tomo_file = it.rlnMicrographName - base_name = os.path.basename(tomo_file) + base_name = os.path.basename(tomo_file) deconv_tomo_name = '{}/{}'.format(deconv_folder,base_name) - + deconv_one(it.rlnMicrographName,deconv_tomo_name,defocus=it.rlnDefocus/10000.0, pixel_size=it.rlnPixelSize,snrfalloff=it.rlnSnrFalloff, deconvstrength=it.rlnDeconvStrength,highpassnyquist=highpassnyquist,tile=tile,ncpu=ncpu) md._setItemValue(it,Label('rlnDeconvTomoName'),deconv_tomo_name) md.write(star_file) logging.info('\n######Isonet done ctf deconvolve######\n') - + except Exception: error_text = traceback.format_exc() f =open('log.txt','a+') @@ -167,9 +166,9 @@ def deconv(self, star_file: str, f.close() logging.error(error_text) - def make_mask(self,star_file, - mask_folder: str = 'mask', - patch_size: int=4, + def make_mask(self,star_file, + mask_folder: str = 'mask', + patch_size: int=4, density_percentage: int=None, std_percentage: int=None, use_deconv_tomo:bool=True, @@ -180,16 +179,16 @@ def make_mask(self,star_file, isonet.py make_mask star_file [--mask_folder] [--patch_size] [--density_percentage] [--std_percentage] [--use_deconv_tomo] [--tomo_idx] :param star_file: path to the tomogram or tomogram folder :param mask_folder: path and name of the mask to save as - :param patch_size: (4) The size of the box from which the max-filter and std-filter are calculated. - :param density_percentage: (50) The approximate percentage of pixels to keep based on their local pixel density. - If this value is not set, the program will look for the parameter in the star file. + :param patch_size: (4) The size of the box from which the max-filter and std-filter are calculated. + :param density_percentage: (50) The approximate percentage of pixels to keep based on their local pixel density. + If this value is not set, the program will look for the parameter in the star file. If this value is not set and not found in star file, the default value 50 will be used. - :param std_percentage: (50) The approximate percentage of pixels to keel based on their local standard deviation. - If this value is not set, the program will look for the parameter in the star file. + :param std_percentage: (50) The approximate percentage of pixels to keel based on their local standard deviation. + If this value is not set, the program will look for the parameter in the star file. If this value is not set and not found in star file, the default value 50 will be used. - :param use_deconv_tomo: (True) If CTF deconvolved tomogram is found in tomogram.star, use that tomogram instead. - :param z_crop: If exclude the top and bottom regions of tomograms along z axis. For example, "--z_crop 0.2" will mask out the top 20% and bottom 20% region along z axis. - :param tomo_idx: (None) If this value is set, process only the tomograms listed in this index. e.g. 1,2,4 or 5-10,15,16 + :param use_deconv_tomo: (True) If CTF deconvolved tomogram is found in tomogram.star, use that tomogram instead. + :param z_crop: If exclude the top and bottom regions of tomograms along z axis. For example, "--z_crop 0.2" will mask out the top 20% and bottom 20% region along z axis. + :param tomo_idx: (None) If this value is set, process only the tomograms listed in this index. e.g. 1,2,4 or 5-10,15,16 """ from IsoNet.bin.make_mask import make_mask logging.basicConfig(format='%(asctime)s, %(levelname)-8s %(message)s', @@ -201,7 +200,7 @@ def make_mask(self,star_file, # write star percentile threshold md = MetaData() md.read(star_file) - if not 'rlnMaskDensityPercentage' in md.getLabels(): + if not 'rlnMaskDensityPercentage' in md.getLabels(): md.addLabels('rlnMaskDensityPercentage','rlnMaskStdPercentage','rlnMaskName') for it in md: md._setItemValue(it,Label('rlnMaskDensityPercentage'),50) @@ -230,7 +229,7 @@ def make_mask(self,star_file, percentile=it.rlnMaskDensityPercentage, threshold=it.rlnMaskStdPercentage, surface = z_crop) - + md._setItemValue(it,Label('rlnMaskName'),mask_out_name) md.write(star_file) logging.info('\n######Isonet done making mask######\n') @@ -259,11 +258,11 @@ def extract(self, :param subtomo_star: (subtomo.star) star file for output subtomograms. :param cube_size: (64) Size of cubes for training, should be divisible by 8, eg. 32, 64. The actual sizes of extracted subtomograms are 1.5 times of this value. :param log_level: ("info") level of the output, either "info" or "debug" - :param use_deconv_tomo: (True) If CTF deconvolved tomogram is found in tomogram.star, use that tomogram instead. + :param use_deconv_tomo: (True) If CTF deconvolved tomogram is found in tomogram.star, use that tomogram instead. """ d = locals() d_args = Arg(d) - + if d_args.log_level == "debug": logging.basicConfig(format='%(asctime)s, %(levelname)-8s [%(filename)s:%(lineno)d] %(message)s' ,datefmt="%H:%M:%S",level=logging.DEBUG,handlers=[logging.StreamHandler(sys.stdout)]) @@ -331,7 +330,7 @@ def refine(self, :param iterations: (30) Number of training iterations. :param data_dir: (data) Temperary folder to save the generated data used for training. :param log_level: (info) debug level, could be 'info' or 'debug' - :param continue_from: (None) A Json file to continue from. That json file is generated at each iteration of refine. + :param continue_from: (None) A Json file to continue from. That json file is generated at each iteration of refine. :param result_dir: ('results') The name of directory to save refined neural network models and subtomograms :param preprocessing_ncpus: (16) Number of cpu for preprocessing. @@ -343,9 +342,9 @@ def refine(self, ************************Denoise settings************************ - :param noise_level: (0.05,0.1,0.15,0.2) Level of noise STD(added noise)/STD(data) after the iteration defined in noise_start_iter. + :param noise_level: (0.05,0.1,0.15,0.2) Level of noise STD(added noise)/STD(data) after the iteration defined in noise_start_iter. :param noise_start_iter: (11,16,21,26) Iteration that start to add noise of corresponding noise level. - :param noise_mode: (None) Filter names when generating noise volumes, can be 'ramp', 'hamming' and 'noFilter' + :param noise_mode: (None) Filter names when generating noise volumes, can be 'ramp', 'hamming' and 'noFilter' ************************Network settings************************ :param drop_out: (0.3) Drop out rate to reduce overfitting. @@ -379,15 +378,15 @@ def predict(self, star_file: str, model: str, output_dir: str='./corrected_tomos :param batch_size: The batch size of the cubes grouped into for network predicting :param normalize_percentile: (True) if normalize the tomograms by percentile. Should be the same with that in refine parameter. :param log_level: ("debug") level of message to be displayed, could be 'info' or 'debug' - :param tomo_idx: (None) If this value is set, process only the tomograms listed in this index. e.g. 1,2,4 or 5-10,15,16 - :param use_deconv_tomo: (True) If CTF deconvolved tomogram is found in tomogram.star, use that tomogram instead. + :param tomo_idx: (None) If this value is set, process only the tomograms listed in this index. e.g. 1,2,4 or 5-10,15,16 + :param use_deconv_tomo: (True) If CTF deconvolved tomogram is found in tomogram.star, use that tomogram instead. :raises: AttributeError, KeyError """ d = locals() d_args = Arg(d) from IsoNet.bin.predict import predict - - + + if d_args.log_level == "debug": logging.basicConfig(format='%(asctime)s, %(levelname)-8s %(message)s', datefmt="%m-%d %H:%M:%S",level=logging.DEBUG,handlers=[logging.StreamHandler(sys.stdout)]) @@ -402,15 +401,15 @@ def predict(self, star_file: str, model: str, output_dir: str='./corrected_tomos f.write(error_text) f.close() logging.error(error_text) - + def check(self): from IsoNet.bin.predict import predict from IsoNet.bin.refine import run import skimage - import PyQt5 + import PyQt5 import tqdm print('IsoNet --version 0.1 installed') - + def gui(self): import IsoNet.gui.Isonet_star_app as app app.main() @@ -425,7 +424,7 @@ def pool_process(p_func,chunks_list,ncpu): # results = p.map(partial_func,chunks_gpu_num_list,chunksize=1) results = list(p.map(p_func,chunks_list)) # return results - + if __name__ == "__main__": core.Display = Display # logging.basicConfig(format='%(asctime)s, %(levelname)-8s %(message)s',datefmt="%m-%d %H:%M:%S",level=logging.INFO) diff --git a/bin/predict.py b/bin/predict.py index 14a7d39..93399dc 100644 --- a/bin/predict.py +++ b/bin/predict.py @@ -24,7 +24,7 @@ def predict(args): logger = logging.getLogger('predict') if args.log_level == "debug": - logging.basicConfig(format='%(asctime)s, %(levelname)-8s [%(filename)s:%(lineno)d] %(message)s', + logging.basicConfig(format='%(asctime)s, %(levelname)-8s [%(filename)s:%(lineno)d] %(message)s', datefmt="%H:%M:%S",level=logging.DEBUG,handlers=[logging.StreamHandler(sys.stdout)]) else: logging.basicConfig(format='%(asctime)s, %(levelname)-8s %(message)s', @@ -37,7 +37,7 @@ def predict(args): args.batch_size = max(4, 2 * args.ngpus) #print('batch_size',args.batch_size) os.environ["CUDA_DEVICE_ORDER"]="PCI_BUS_ID" - os.environ["CUDA_VISIBLE_DEVICES"]=args.gpuID + os.environ["CUDA_VISIBLE_DEVICES"]=args.gpuID #check gpu settings from IsoNet.bin.refine import check_gpu check_gpu(args) @@ -59,7 +59,7 @@ def predict(args): # write star percentile threshold md = MetaData() md.read(args.star_file) - if not 'rlnCorrectedTomoName' in md.getLabels(): + if not 'rlnCorrectedTomoName' in md.getLabels(): md.addLabels('rlnCorrectedTomoName') for it in md: md._setItemValue(it,Label('rlnCorrectedTomoName'),None) @@ -76,7 +76,7 @@ def predict(args): predict_one(args,tomo_file,model,output_file=tomo_out_name) md._setItemValue(it,Label('rlnCorrectedTomoName'),tomo_out_name) md.write(args.star_file) - + def predict_one(args,one_tomo,model,output_file=None): #predict one tomogram in mrc format INPUT: mrc_file string OUTPUT: output_file(str) or _corrected.mrc import logging @@ -92,6 +92,7 @@ def predict_one(args,one_tomo,model,output_file=None): with mrcfile.open(one_tomo) as mrcData: real_data = mrcData.data.astype(np.float32)*-1 + voxelsize = mrcData.voxel_size real_data = normalize(real_data,percentile=args.normalize_percentile) data=np.expand_dims(real_data,axis=-1) reform_ins = reform3D(data) @@ -100,7 +101,7 @@ def predict_one(args,one_tomo,model,output_file=None): #imposing wedge to every cubes #data=wedge_imposing(data) - N = args.batch_size * args.ngpus * 4 # 8*4*8 + N = args.batch_size * args.ngpus * 4 # 8*4*8 num_patches = data.shape[0] if num_patches%N == 0: append_number = 0 @@ -122,5 +123,6 @@ def predict_one(args,one_tomo,model,output_file=None): outData = normalize(outData,percentile=args.normalize_percentile) with mrcfile.new(output_file, overwrite=True) as output_mrc: output_mrc.set_data(-outData) + output_mrc.voxel_size = voxelsize logging.info('Done predicting') - # predict(args.model,args.weight,args.mrc_file,args.output_file, cubesize=args.cubesize, cropsize=args.cropsize, batch_size=args.batch_size, gpuID=args.gpuID, if_percentile=if_percentile) \ No newline at end of file + # predict(args.model,args.weight,args.mrc_file,args.output_file, cubesize=args.cubesize, cropsize=args.cropsize, batch_size=args.batch_size, gpuID=args.gpuID, if_percentile=if_percentile) diff --git a/util/deconvolution.py b/util/deconvolution.py index 560c024..d329a38 100644 --- a/util/deconvolution.py +++ b/util/deconvolution.py @@ -44,6 +44,7 @@ def wiener1d(angpix, defocus, snrfalloff, deconvstrength, highpassnyquist, phase def tom_deconv_tomo(vol_file, angpix, defocus, snrfalloff, deconvstrength, highpassnyquist, phaseflipped, phaseshift): with mrcfile.open(vol_file) as f: vol = f.data + voxelsize = f.voxel_size data = np.arange(0,1+1/2047.,1/2047.) highpass = np.minimum(np.ones(data.shape[0]), data/highpassnyquist) * np.pi; highpass = 1-np.cos(highpass); @@ -65,15 +66,15 @@ def tom_deconv_tomo(vol_file, angpix, defocus, snrfalloff, deconvstrength, highp s1 = - int(np.shape(vol)[1] / 2) f1 = s1 + np.shape(vol)[1] - 1 - m1 = np.arange(s1,f1+1) + m1 = np.arange(s1,f1+1) s2 = - int(np.shape(vol)[0] / 2) f2 = s2 + np.shape(vol)[0] - 1 - m2 = np.arange(s2,f2+1) + m2 = np.arange(s2,f2+1) s3 = - int(np.shape(vol)[2] / 2) f3 = s3 + np.shape(vol)[2] - 1 - m3 = np.arange(s3,f3+1) + m3 = np.arange(s3,f3+1) #s3 = -floor(size(vol,3)/2); #f3 = s3 + size(vol,3) - 1; @@ -97,6 +98,7 @@ def tom_deconv_tomo(vol_file, angpix, defocus, snrfalloff, deconvstrength, highp deconv = deconv/np.std(deconv) * np.std(vol) + np.average(vol) with mrcfile.new(os.path.splitext(vol_file)[0]+'_deconv.mrc',overwrite=True) as n: n.set_data(deconv.astype(np.float32)) #.astype(type(vol[0,0,0])) + n.voxel_size = voxelsize #return real(ifftn(fftn(single(vol)).*ramp)); return os.path.splitext(vol_file)[0]+'_deconv.mrc' @@ -171,7 +173,7 @@ def restore(self,new_file_list): new = np.zeros((self._N[0]*cubesize,self._N[1]*cubesize,self._N[2]*cubesize),dtype = np.float32) start=int((cropsize-cubesize)/2) end=int((cropsize+cubesize)/2) - + for i in range(self._N[0]): for j in range(self._N[1]): for k in range(self._N[2]): @@ -194,10 +196,10 @@ def restore(self,new_file_list): one_chuck = f.data # print(one_chuck[overlap_len[0]//2:-(overlap_len[0]//2),overlap_len[1]//2:-(overlap_len[1]//2),overlap_len[2]//2:-(overlap_len[2]//2)].shape) # print(one_chuck.shape) - new_vol[i[0]+overlap_len[0]//2:i[1]-overlap_len[0]//2,j[0]+overlap_len[1]//2:j[1]-overlap_len[1]//2, + new_vol[i[0]+overlap_len[0]//2:i[1]-overlap_len[0]//2,j[0]+overlap_len[1]//2:j[1]-overlap_len[1]//2, k[0]+overlap_len[2]//2:k[1]-overlap_len[2]//2] = one_chuck[overlap_len[0]//2:-(overlap_len[0]//2),overlap_len[1]//2:-(overlap_len[1]//2),overlap_len[2]//2:-(overlap_len[2]//2)] - - + + # return np.multiply(new_vol,1/factor_vol) return new_vol[overlap_len[0]//2:-(overlap_len[0]//2), overlap_len[1]//2:-(overlap_len[1]//2), @@ -214,20 +216,20 @@ def deconv_one(tomo, out_tomo,defocus=1.0, pixel_size=1.0,snrfalloff=1.0, deconv :param ngpu: (4) number of avaliable gpu cards :param gpu_memory: (10) memory of each gpu :param pixel_size: (10) pixel size in anstroms - :param: snrfalloff: (1.0) The larger this values, more high frequency informetion are filtered out. - :param deconvstrength: (1.0) + :param: snrfalloff: (1.0) The larger this values, more high frequency informetion are filtered out. + :param deconvstrength: (1.0) """ import mrcfile from multiprocessing import Pool from functools import partial from IsoNet.util.deconvolution import tom_deconv_tomo,Chunks import shutil - import time + import time t1 = time.time() if os.path.isdir('./deconv_temp'): shutil.rmtree('./deconv_temp') os.mkdir('./deconv_temp') - + root_name = os.path.splitext(os.path.basename(tomo))[0] logging.info('deconv: {}| pixel: {}| defocus: {}| snrfalloff:{}| deconvstrength:{}'.format(tomo, pixel_size, defocus ,snrfalloff, deconvstrength)) @@ -237,19 +239,23 @@ def deconv_one(tomo, out_tomo,defocus=1.0, pixel_size=1.0,snrfalloff=1.0, deconv # print(chunks_list) chunks_deconv_list = [] with Pool(ncpu) as p: - partial_func = partial(tom_deconv_tomo,angpix=pixel_size, defocus=defocus, snrfalloff=snrfalloff, + partial_func = partial(tom_deconv_tomo,angpix=pixel_size, defocus=defocus, snrfalloff=snrfalloff, deconvstrength=deconvstrength, highpassnyquist=highpassnyquist, phaseflipped=False, phaseshift=0 ) # results = p.map(partial_func,chunks_gpu_num_list,chunksize=1) chunks_deconv_list = list(p.map(partial_func,chunks_list)) # pool_process(partial_func,chunks_list_single_pool,ncpu) # chunks_deconv_list += results vol_restored = c.restore(chunks_deconv_list) + with mrcfile.open(tomo) as mrc: + voxelsize = mrc.voxel_size + with mrcfile.new(out_tomo, overwrite=True) as mrc: mrc.set_data(vol_restored) + mrc.voxel_size = voxelsize shutil.rmtree('./deconv_temp') t2 = time.time() logging.info('time consumed: {:10.4f} s'.format(t2-t1)) - + @@ -261,16 +267,17 @@ def deconv_gpu(tomo, defocus: float=1.0, pixel_size: float=1.0,snrfalloff: float :param tomo: tomogram file :param defocus: (1) defocus in um :param pixel_size: (10) pixel size in anstroms - :param: snrfalloff: (1.0) The larger this values, more high frequency informetion are filtered out. - :param deconvstrength: (1.0) + :param: snrfalloff: (1.0) The larger this values, more high frequency informetion are filtered out. + :param deconvstrength: (1.0) """ import mrcfile from multiprocessing import Pool from functools import partial from IsoNet.util.deconv_gpu import Chunks,tom_deconv_tomo import sys - # from IsoNet.util.deconvolution import + # from IsoNet.util.deconvolution import with mrcfile.open(tomo) as mrc: + voxelsize = mrc.voxel_size vol = mrc.data outname = os.path.splitext(os.path.basename(tomo))[0] +'-deconv.rec' @@ -281,7 +288,7 @@ def deconv_gpu(tomo, defocus: float=1.0, pixel_size: float=1.0,snrfalloff: float chunks_gpu_num_list = [[array,j%num_gpu] for j,array in enumerate(chunks_list)] print('chunks_list',chunks_list.__sizeof__()) with Pool(ncpu) as p: - partial_func = partial(tom_deconv_tomo,angpix=pixel_size, defocus=defocus, snrfalloff=snrfalloff, + partial_func = partial(tom_deconv_tomo,angpix=pixel_size, defocus=defocus, snrfalloff=snrfalloff, deconvstrength=deconvstrength, highpassnyquist=0.1, phaseflipped=False, phaseshift=0 ) # chunks_deconv_list = list(p.map(partial_func,chunks_gpu_num_list,chunksize=1)) # results = p.map(partial_func,chunks_list) @@ -309,5 +316,5 @@ def deconv_gpu(tomo, defocus: float=1.0, pixel_size: float=1.0,snrfalloff: float parser.add_argument("--ncpu",type=int,default=8) args = parser.parse_args() start = time.time() - + deconv_one(args.mrcFile, args.outFile,defocus=args.defocus/10000.0, pixel_size=args.pixsize,snrfalloff=args.snrfalloff, deconvstrength=args.deconvstrength,tile=args.tile,ncpu=args.ncpu)