cryoSPARC调用motioncorr2程序运行报错
报错信息:
Traceback (most recent call last):
File "cryosparc_master/cryosparc_compute/run.py", line 116, in cryosparc_master.cryosparc_compute.run.main
File "/root/cryosparcuser/cryosparc/cryosparc_worker/cryosparc_compute/jobs/motioncorrection/run_motioncor2.py", line 398, in run_motioncor2_wrapper
global_shifts = parse_global_alignment(full_frame_log_abs)
File "/root/cryosparcuser/cryosparc/cryosparc_worker/cryosparc_compute/jobs/motioncorrection/run_motioncor2.py", line 22, in parse_global_alignment
with open(full_frame_alignment_log) as logfile:
FileNotFoundError: [Errno 2] No such file or directory: '/data6/x/CS-x/J72/motioncor2_logs/00-Patch-Full.log'
原因
worker安装路径下/cryosparc_compute/jobs/motioncorrection/run_motioncor2.py
代码与motioncorr2_1.6.3版本调用命令不同
解决方案:
替换run_motioncor2.py代码内容
## ---------------------------------------------------------------------------
## Copyright (c) 2024 Structura Biotechnology Inc. All rights reserved.
## Do not reproduce or redistribute, in whole or in part.
## Use of this code is permitted only under licence from Structura.
## Contact us at info@structura.bio.
## ---------------------------------------------------------------------------junquan
from builtins import str
from builtins import range
import os
import time
import numpy as n
from .. import runcommon as rc
import matplotlib.pyplot as plt
import matplotlib
from ...blobio import mrc
from ...blobio import tiff
def parse_global_alignment(full_frame_alignment_log):
shifts = []
with open(full_frame_alignment_log) as logfile:
for line in logfile:
stripped_line = line.strip()
if '#' not in stripped_line and len(stripped_line) > 0:
values = stripped_line.split()
shifts.append([float(values[1]), float(values[2])])
global_shifts = n.array([shifts])
return global_shifts
def parse_patch_alignment(patch_patch_alignment_log):
shifts = []
patches = []
with open(patch_patch_alignment_log) as logfile:
for line in logfile:
stripped_line = line.strip()
if '#' in stripped_line or len(stripped_line) == 0: continue
values = stripped_line.split()
if int(values[0]) == 1:
# e.g. 1 740.91 766.51 0.00 0.00 0.00 0.00
# this means we've reached a new patch
shifts.append([])
# record the coordinates of the patch
patches.append([float(values[1]), float(values[2])])
shifts[-1].append([float(values[-2]), float(values[-1])])
patch_shifts = n.array(shifts)
patch_coords = n.array(patches)
return patch_shifts, patch_coords
def find_extension(abs_path):
ext = str(os.path.splitext(abs_path)[1])
if ext in ['.mrc', '.mrcs', '.stk']:
return 'mrc'
elif ext in ['.tiff', '.tif']:
return 'tif'
else:
assert False, "Cant read extension %s for input file path" % (ext)
def make_command(
exec_path,
abs_input_path,
abs_output_path,
abs_logfile_path,
gpu_resources,
gpu_mem_usage,
patch_x,
patch_y,
patch_overlap,
abs_gain_ref_path,
gain_rotate_num,
gain_flip_num,
abs_dark_ref_path,
abs_defect_correct_path,
alignment_config,
iter,
tol,
throw,
trunc,
align_to_frame_num,
dose_weight,
psize_A,
accel_kv,
frame_dose,
tilt_start,
tilt_step,
selected_sum,
sum_range_a,
sum_range_b,
major_scale,
minor_scale,
distortion_angle,
fourier_cropping,
low_snr,
b_factor_global,
b_factor_local,
group_num,
frame_motion_blurring,
phase_only,
extra
):
path_to_motioncor2_executable = exec_path
command = [
path_to_motioncor2_executable,
]
#find extension of input movie
input_type = find_extension(abs_input_path)
if input_type == "tif": command.extend(['-InTiff', abs_input_path])
elif input_type == "mrc": command.extend(['-InMrc', abs_input_path])
if abs_output_path is not None: command.extend(['-OutMrc', abs_output_path])
if abs_gain_ref_path is not None: command.extend(['-Gain', abs_gain_ref_path])
if gain_rotate_num != 0: command.extend(['-RotGain', gain_rotate_num])
if gain_flip_num != 0: command.extend(['-FlipGain', gain_flip_num])
if abs_dark_ref_path is not None:
dark_ext = str(os.path.splitext(abs_dark_ref_path)[1])
assert dark_ext in ['.mrc', '.mrcs', '.stk'], "Dark Reference file must be a .mrc file"
command.extend(['-Dark', abs_dark_ref_path])
if abs_defect_correct_path is not None:
defect_ext = str(os.path.splitext(abs_defect_correct_path)[1])
assert defect_ext in ['.txt'], "Defect file must be a .txt file"
command.extend(['-DefectFile', abs_defect_correct_path])
if patch_x is not None and patch_y is not None:
if patch_overlap is None:
command.extend(['-Patch', patch_x, patch_y])
else:
command.extend(['-Patch', patch_x, patch_y, patch_overlap])
if alignment_config: command.extend(['-Iter', iter, '-Tol', tol])
if throw is not None: command.extend(['-Throw', throw])
if trunc is not None: command.extend(['-Trunc', trunc])
if align_to_frame_num is not None: command.extend(['-FmRef', align_to_frame_num])
if dose_weight: command.extend(['-Kv', accel_kv, '-PixSize', psize_A, '-FmDose', frame_dose])
if (tilt_start is not None and tilt_step is not None): command.extend(['-Tilt', tilt_start, tilt_step])
if selected_sum: command.extend(['-SumRange', sum_range_a, sum_range_b])
if major_scale is not None and minor_scale is not None and distortion_angle is not None: command.extend(['-Mag', major_scale, minor_scale, distortion_angle])
if fourier_cropping is not None: command.extend(['-FtBin', fourier_cropping])
if low_snr:
if (b_factor_global is not None) & (b_factor_local is not None): command.extend(['-Bft', b_factor_global, b_factor_local])
if group_num is not None: command.extend (['-Group', group_num])
if frame_motion_blurring: command.extend(['-InFmMotion 1'])
if phase_only: command.extend(['-PhaseOnly'])
if gpu_resources != 0:
command.extend(['-Gpu', gpu_resources])
command.extend(['-GpuMemUsage', gpu_mem_usage])
command.extend(['-LogDir', abs_logfile_path]) #---------------------------------junquan-----------------------
if extra is not None: command.extend([extra])
command = [str(c) for c in command] # convert everything to string
return(command)
def run_motioncor2_wrapper(job):
""" function that runs the actual job. """
params = rc.com.get_merged_params(job)
proj_dir_abs = rc.get_project_dir_abs()
job_dir_rel = job['job_dir']
job_dir_abs = os.path.join(proj_dir_abs, job_dir_rel)
#puid = job['project_uid']
#juid = job['uid']
res_alloc = job['resources_allocated']
hostname = res_alloc['hostname']
print("Running job on hostname %s", hostname)
print("Allocated Resources : ", res_alloc)
cuda_devs = res_alloc['slots']['GPU']
cuda_dev = ("".join(str(dev) + " " for dev in cuda_devs)).strip() #converts [0,1,2] into '0 1 2'
#terms
rc.log("CryoSPARC does not distribute MotionCor2 binaries. Please ensure you have your own copy of MotionCor2 installed under the terms of the MotionCor2 Non-Commercial Software License Agreement available at: https://docs.google.com/forms/d/e/1FAIpQLSfAQm5MA81qTx90W9JL6ClzSrM77tytsvyyHh1ZZWrFByhmfQ/formResponse. For-profit users must contact David Agard for licensing information prior to download. Structura Biotechnology Inc. makes no warranty regarding MotionCor2.")
movies_dset = rc.load_input_group(input_group_name='movies')
movies_dset.add_fields(rc.com.get_dataset_columns('exposure.micrograph_blob', 'micrograph_blob_non_dw'))
movies_dset.add_fields(rc.com.get_dataset_columns('exposure.motion', 'rigid_motion'))
movies_dset.add_fields(rc.com.get_dataset_columns('exposure.motion', 'local_motion'))
#if the user selected the dose weighting option (which applies to all micrographs outputted)
do_dose_weighting = params.get('dose_weighting', False)
if do_dose_weighting:
movies_dset.add_fields(rc.com.get_dataset_columns('exposure.micrograph_blob', 'micrograph_blob'))
os.makedirs(os.path.join(proj_dir_abs, job_dir_rel, 'motioncorrected'))
rc.log("Created dirctory %s" %(os.path.join(proj_dir_abs, job_dir_rel, 'motioncorrected')))
os.makedirs(os.path.join(proj_dir_abs, job_dir_rel, 'motioncor2_logs'))
rc.log("Created directory %s"% (os.path.join(proj_dir_abs, job_dir_rel, 'motioncor2_logs')))
#TODO: Find a better way to store the path to motioncor2 executable (maybe in config.sh?)
motioncor2_exec_path = params.get('exec_path', None)
assert motioncor2_exec_path is not None, "Must provide path to motioncor2 executable."
movies = movies_dset.rows()
keep_idxs = []
#TODO: BATCH PROCESSING MODE
# SINGLE PROCESSING MODE
num_mov_processed = 0
checkpoint_freq = 5
num_mov_total = len(movies[:params.get('random_num')])
fulltic = time.time()
rc.log("Starting to process exposures...")
for movie in movies[:params.get('random_num')]:
if num_mov_processed % checkpoint_freq == 0:
rc.log_checkpoint()
rc.log("--------------------------------------------------------------")
rc.log("Processed %d of %d movies in %.2fs " % (num_mov_processed, num_mov_total, time.time()-fulltic))
keep_idxs.append(movie.idx)
movietic = time.time()
tic = time.time()
#find the full path to the input movie
mov_path_rel = movie['movie_blob/path']
mov_path_abs = os.path.join(proj_dir_abs, mov_path_rel)
rc.log('Raw movie filepath located at: %s - creating MotionCor2 command string...' % (mov_path_rel))
input_type = find_extension(mov_path_abs)
if input_type == "mrc":
mov_raw = mrc.read_mrc(mov_path_abs)
elif input_type == "tif":
mov_raw = tiff.read_tiff(mov_path_abs)
#full path to output micrograph
basename = rc.com.get_movie_basename(movie['movie_blob/path'])
output_path_rel = os.path.join(job_dir_rel, 'motioncorrected', basename + '_motioncor2_aligned.mrc' )
output_path_abs = os.path.join(proj_dir_abs, output_path_rel)
#full path to dose weighted micrograph
if do_dose_weighting:
dw_output_path_rel = os.path.join(job_dir_rel, 'motioncorrected', basename + '_motioncor2_aligned_DW.mrc' )
dw_output_path_abs = os.path.join(proj_dir_abs, dw_output_path_rel)
#motioncor2 log prefix file
motioncor2_log_file_prefix = os.path.join(proj_dir_abs, job_dir_rel, 'motioncor2_logs') #, str(movie.idx)
psize_raw = movie['movie_blob/psize_A']
psize_final = psize_raw
fourier_cropping = params.get('fourier_cropping', None)
if fourier_cropping is not None:
psize_final = psize_raw * fourier_cropping
#obtain movie's microscope params (its a requirement that these are provided)
if not params['override_dose_weighting']:
accel_kv = movie['mscope_params/accel_kv']
movie_shape = tuple(movie['movie_blob/shape'])
frame_dose = movie['mscope_params/total_dose_e_per_A2'] / movie_shape[0] #divide frame dose by total number of frames
else:
accel_kv = params['accel_kv']
frame_dose = params['fm_dose_e_per_A2']
#get gain reference details
gain_path_abs = None; flip_x = False; flip_y = False; rotate_num = 0; flip_num = 0
if not movie['movie_blob/is_gain_corrected']:
gain_path_rel = movie['gain_ref_blob/path']
gain_path_abs = os.path.join(proj_dir_abs, gain_path_rel)
flip_x = movie['gain_ref_blob/flip_x']
flip_y = movie['gain_ref_blob/flip_y']
if flip_x:
flip_num = 2
elif flip_y:
flip_num = 1
gain_rotate_num = movie['gain_ref_blob/rotate_num']
if (gain_rotate_num >= 1) & (gain_rotate_num <= 91):
rotate_num = 1
elif (gain_rotate_num > 91) & (gain_rotate_num <= 181):
rotate_num = 2
elif (gain_rotate_num > 181) & (gain_rotate_num <= 271):
rotate_num = 3
command = make_command( exec_path = motioncor2_exec_path,
abs_input_path = mov_path_abs,
abs_output_path = output_path_abs,
abs_logfile_path = motioncor2_log_file_prefix,
gpu_resources = cuda_dev,
gpu_mem_usage = params.get('gpu_mem_usage', 0.5),
patch_x = params.get('patch_x', 5.0),
patch_y = params.get('patch_y', 5.0),
patch_overlap = params.get('patch_overlap', None),
abs_gain_ref_path = gain_path_abs,
gain_rotate_num = rotate_num,
gain_flip_num = flip_num,
abs_dark_ref_path = params.get('dark_ref_path', None),
abs_defect_correct_path = params.get('defect_path', None),
alignment_config = params.get('alignment_config', False),
iter = params.get('max_iterations', None),
tol = params.get('alignment_tolerance', None),
throw = params.get('throw', None),
trunc = params.get('trunc', None),
align_to_frame_num = params.get('frame_index', None),
dose_weight = params.get('dose_weighting', False),
psize_A = psize_raw,
accel_kv = accel_kv,
frame_dose = frame_dose,
tilt_start = params.get('tilt_start', None),
tilt_step = params.get('tilt_step', None),
selected_sum = params.get('selected_sum', False),
sum_range_a = params.get('sum_range_a', None),
sum_range_b = params.get('sum_range_b', None),
major_scale = params.get('major_scale', None),
minor_scale = params.get('minor_scale', None),
distortion_angle = params.get('distortion_angle', None),
fourier_cropping = fourier_cropping,
low_snr = params.get('low_snr', False),
b_factor_global = params.get('b_factor_global', None),
b_factor_local = params.get('b_factor_local', None),
group_num = params.get('group_num', 1),
frame_motion_blurring = params.get('frame_motion_blurring', False),
phase_only = params.get('phase_only', False),
extra = params.get('extra', None),
)
rc.log('Finished creating MotionCor2 command string in %.2fs' % (time.time() - tic))
rc.log('Starting MotionCor2 process...')
tic = time.time()
try:
rc.log('Running MotionCor2 command:', (' ').join(command))
joblog = job_dir_abs + '/job.log'
process = rc.start_external_process(command, log_file_path=joblog)
rc.log('Running process %i' % process.pid)
process.communicate()
rc.remove_subprocess(process.pid)
assert os.path.exists(output_path_abs), 'motioncor2 failed to produce output file %s' % output_path_abs
except Exception as e:
rc.log("ERROR", e)
rc.log('Finished MotionCor2 process in %.2fs' % (time.time() - tic))
tic = time.time()
with open(output_path_abs) as mrc_file:
output_mic_header = mrc.read_mrc_header(mrc_file, debug_fname=output_path_abs)
movie['micrograph_blob_non_dw/path'] = output_path_rel
movie['micrograph_blob_non_dw/idx'] = 0
movie['micrograph_blob_non_dw/shape'] = [output_mic_header['ny'], output_mic_header['nx']]
movie['micrograph_blob_non_dw/psize_A'] = psize_final
movie['micrograph_blob_non_dw/format'] = 'MRC/2'
movie['micrograph_blob_non_dw/is_background_subtracted'] = False
if do_dose_weighting:
with open(dw_output_path_abs) as mrc_file:
dw_output_mic_header = mrc.read_mrc_header(mrc_file, debug_fname=dw_output_path_abs)
movie['micrograph_blob/path'] = dw_output_path_rel
movie['micrograph_blob/idx'] = 0
movie['micrograph_blob/shape'] = [dw_output_mic_header['ny'], dw_output_mic_header['nx']]
movie['micrograph_blob/psize_A'] = psize_final
movie['micrograph_blob/format'] = 'MRC/2'
movie['micrograph_blob/is_background_subtracted'] = False
frame_start = params.get('throw', 0)
frame_end = params.get('trunc', None)
zero_shift_frame = params.get('align_to_frame_num', 0)
movie['rigid_motion/type'] = 'rigid'
movie['rigid_motion/idx'] = 0
movie['rigid_motion/frame_start'] = 0 if frame_start is None else frame_start
movie['rigid_motion/frame_end'] = len(mov_raw) if frame_end is None else frame_end
movie['rigid_motion/zero_shift_frame'] = 0 if zero_shift_frame is None else zero_shift_frame # BEFORE frame start is applied
rigid_motion_traj_rel = os.path.join(job_dir_rel, 'motioncorrected', basename + '_rigid_motion_traj.npy' )
rigid_motion_traj_abs = os.path.join(proj_dir_abs, rigid_motion_traj_rel)
full_frame_log_abs = motioncor2_log_file_prefix+'/'+basename+'-Patch-Full.log'
global_shifts = parse_global_alignment(full_frame_log_abs)
rc.com.write_array_to_numpy_file(global_shifts, rigid_motion_traj_abs)
movie['rigid_motion/path'] = rigid_motion_traj_rel
movie['local_motion/type'] = 'patch'
movie['local_motion/idx'] = 0
movie['local_motion/frame_start'] = 0 if frame_start is None else frame_start
movie['local_motion/frame_end'] = len(mov_raw) if frame_end is None else frame_end
movie['local_motion/zero_shift_frame'] = 0 if zero_shift_frame is None else zero_shift_frame # BEFORE frame start is applied
local_patch_traj_rel = os.path.join(job_dir_rel, 'motioncorrected', basename + '_local_patch_traj.npy' )
local_patch_traj_abs = os.path.join(proj_dir_abs, local_patch_traj_rel)
patch_patch_log_abs = motioncor2_log_file_prefix+'/'+basename+'-Patch-Patch.log'
patch_shifts, patch_coords = parse_patch_alignment(patch_patch_log_abs)
rc.com.write_array_to_numpy_file(patch_shifts, local_patch_traj_abs)
movie['local_motion/path'] = local_patch_traj_rel
if num_mov_processed < params['num_plots']:
fit_cumulative_patch_shifts = []
for patch in range(0,len(patch_coords)):
fit_cumulative_patch_shifts.append([])
for shift in patch_shifts[patch]:
fit_cumulative_patch_shifts[-1].append([patch_coords[patch][0]+shift[0], patch_coords[patch][1]+shift[1]])
fit_cumulative_patch_shifts = n.array(fit_cumulative_patch_shifts)
matplotlib.rcParams.update({'font.size': 8})
motioncorrected_figures = plt.figure(figsize=(8,4))
plt.subplot(121)
plt.title('Local Motion Correction (Patch 0)')
plt.plot(fit_cumulative_patch_shifts[0][:,0],fit_cumulative_patch_shifts[0][:,1],'-k.', label='Patch Motion Alignment')
plt.axis('equal')
plt.grid(True)
plt.subplot(122)
plt.title('Global Motion Correction')
plt.plot(global_shifts[0][:,0],global_shifts[0][:,1],'-k.', label='Global Motion Alignment')
plt.axis('equal')
plt.grid(True)
plt.tight_layout()
rc.log_plot(motioncorrected_figures, 'Local and Global Motion for %s' % basename)
if num_mov_processed < 1:
rc.set_tile_image('img%d'%num_mov_processed, motioncorrected_figures, 1, 2)
if num_mov_processed < 1:
rc.set_output_group_image('micrographs', motioncorrected_figures)
rc.log('Finished plotting in %.2fs' % (time.time() - tic))
rc.log('================ Completed motion correction of exposure in %.2fs' % (time.time() - movietic))
num_mov_processed += 1
outpath_rel = os.path.join(job_dir_rel, 'movies_motioncor2_aligned.cs')
movies_dset = movies_dset.take(keep_idxs)
movies_dset.save(os.path.join(proj_dir_abs, outpath_rel))
rc.output('micrographs', 'micrograph_blob_non_dw', outpath_rel, 0, len(movies_dset))
if do_dose_weighting:
rc.output('micrographs', 'micrograph_blob', outpath_rel, 0, len(movies_dset))
rc.output('micrographs', 'rigid_motion', outpath_rel, 0, len(movies_dset))
rc.output('micrographs', 'local_motion', outpath_rel, 0, len(movies_dset))
rc.log("Finished processing all exposures!")