Skip to content

Commit

Permalink
fixed seed during SMI, read input .bshape/.echotime files, drop progr…
Browse files Browse the repository at this point in the history
…ess bars, github actions again
  • Loading branch information
Benjamin Aron authored and Benjamin Aron committed Jun 17, 2024
1 parent 56aeb08 commit 2cc4ae2
Show file tree
Hide file tree
Showing 8 changed files with 215 additions and 77 deletions.
29 changes: 14 additions & 15 deletions .github/workflows/build_wheels.yml
Original file line number Diff line number Diff line change
Expand Up @@ -39,8 +39,11 @@ jobs:
fi
done
"
- name: Verify wheels
run: |
ls -al dist/
- name: Upload to PyPI
uses: pypa/gh-action-pypi-publish@master
uses: pypa/gh-action-pypi-publish@release/v1
with:
repository_url: https://upload.pypi.org/legacy/
username: __token__
Expand Down Expand Up @@ -79,8 +82,11 @@ jobs:
delocate-listdeps dist/*.whl # lists library dependencies
delocate-wheel dist/*.whl # copies library dependencies into the wheel
fi
- name: Verify wheels
run: |
ls -al dist/
- name: Upload to PyPI
uses: pypa/gh-action-pypi-publish@master
uses: pypa/gh-action-pypi-publish@release/v1
with:
repository_url: https://upload.pypi.org/legacy/
username: __token__
Expand All @@ -100,24 +106,17 @@ jobs:
- name: Install dependencies
run: |
python -m pip install --upgrade pip setuptools wheel pybind11
- name: Download and install FFTW
- name: Install FFTW
run: |
curl -L -o fftw.zip http://www.fftw.org/fftw-3.3.10-dll64.zip
if (-not (Test-Path fftw.zip)) {
Write-Error "Failed to download fftw.zip"
exit 1
}
Expand-Archive fftw.zip -DestinationPath fftw
if (-not (Test-Path fftw)) {
Write-Error "Failed to unzip fftw.zip"
exit 1
}
copy fftw/*.dll C:\Windows\System32
conda install -c conda-forge fftw
- name: Build wheels
run: |
python setup.py bdist_wheel
- name: Verify wheels
run: |
ls -al dist/
- name: Upload to PyPI
uses: pypa/gh-action-pypi-publish@master
uses: pypa/gh-action-pypi-publish@release/v1
with:
repository_url: https://upload.pypi.org/legacy/
username: __token__
Expand Down
140 changes: 130 additions & 10 deletions lib/designer_func_wrappers.py
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,9 @@
run_normalization:
"""

import time
import shutil

def run_mppca(args_extent, args_phase, args_shrinkage, args_algorithm):
"""
wrapper for complex or magnitude, adatpive or local mppca
Expand Down Expand Up @@ -45,7 +48,13 @@ def run_mppca(args_extent, args_phase, args_shrinkage, args_algorithm):

# note adaptive patch does not include mpcomplex
n_cores = app.ARGS.n_cores
terminal_width = shutil.get_terminal_size().columns
separator = "=" * terminal_width

print("\n" + separator)
print('...denoising...')
start_time = time.time()

if app.ARGS.adaptive_patch:
Signal, Sigma, Nparameters = mp.denoise(
dwi, phase=args_phase, patchtype='nonlocal', shrinkage=args_shrinkage, algorithm=args_algorithm
Expand Down Expand Up @@ -75,6 +84,17 @@ def run_mppca(args_extent, args_phase, args_shrinkage, args_algorithm):
#run.command('dwidenoise -noise fullnoisemap.mif -estimator Exp2 working.mif dwidn.mif')
run.command('mrconvert -force dwidn.mif working.mif', show=False)

# End timer
end_time = time.time()
elapsed_time = end_time - start_time

# Convert elapsed time to hours and minutes
hours, rem = divmod(elapsed_time, 3600)
minutes, seconds = divmod(rem, 60)

print(f"Denoising completed in {int(hours):02} hours, {int(minutes):02} minutes, {int(seconds):02} seconds.")
print(separator + "\n")

def run_patch2self():
"""
wrapper for patch2self
Expand All @@ -84,6 +104,13 @@ def run_patch2self():
from ants import image_read, image_write, from_numpy
import numpy as np

terminal_width = shutil.get_terminal_size().columns
separator = "=" * terminal_width

print("\n" + separator)
print('...denoising...')
start_time = time.time()

run.command('mrconvert -force -export_grad_fsl working.bvec working.bval working.mif working.nii', show=False)
nii = image_read('working.nii')
dwi = nii.numpy()
Expand All @@ -95,6 +122,16 @@ def run_patch2self():
image_write(out, 'working_p2s.nii')
run.command('mrconvert -force -fslgrad working.bvec working.bval working_p2s.nii working.mif', show=False)

# End timer
end_time = time.time()
elapsed_time = end_time - start_time

# Convert elapsed time to hours and minutes
hours, rem = divmod(elapsed_time, 3600)
minutes, seconds = divmod(rem, 60)

print(f"Denoising completed in {int(hours):02} hours, {int(minutes):02} minutes, {int(seconds):02} seconds.")
print(separator + "\n")

def run_degibbs(pf, pe_dir):
"""
Expand All @@ -104,7 +141,7 @@ def run_degibbs(pf, pe_dir):
import lib.rpg as rpg
from mrtrix3 import run, app, MRtrixError
from ants import image_read, image_write, from_numpy
from tqdm import tqdm
# from tqdm import tqdm
# import os

# rpg_path = os.path.join(os.path.dirname(os.path.dirname(os.path.abspath(__file__))),'rpg_cpp')
Expand All @@ -113,7 +150,13 @@ def run_degibbs(pf, pe_dir):
run.command('mrconvert -force -export_grad_fsl working.bvec working.bval working.mif working.nii', show=False)
nii = image_read('working.nii')
dwi = nii.numpy()
print("...RPG degibbsing...")

terminal_width = shutil.get_terminal_size().columns
separator = "=" * terminal_width

print("\n" + separator)
print('...RPG degibbsing...')
start_time = time.time()

n_cores = app.ARGS.n_cores

Expand All @@ -131,13 +174,13 @@ def run_degibbs(pf, pe_dir):
pe_dir = 1

dwi_t = dwi.transpose(3,2,1,0)
progress_bar = tqdm(total=100)
def progress_callback(progress):
progress_bar.n = progress
progress_bar.refresh()
#progress_bar = tqdm(total=100)
#def progress_callback(progress):
# progress_bar.n = progress
# progress_bar.refresh()

dwi_dg_t = rpg.unring(dwi_t, minW=1, maxW=3, nsh=20, pfv=float(pf), pfdimf=pe_dir, phase_flag=False, progress_callback=progress_callback)
progress_bar.close()
dwi_dg_t = rpg.unring(dwi_t, minW=1, maxW=3, nsh=20, pfv=float(pf), pfdimf=pe_dir, phase_flag=False)
#progress_bar.close()
dwi_dg = dwi_dg_t[0].copy().transpose(3,2,1,0)

out = from_numpy(
Expand All @@ -147,14 +190,29 @@ def progress_callback(progress):
#convert gibbs corrected nii to .mif
run.command('mrconvert -force -fslgrad working.bvec working.bval working_rpg.nii working.mif', show=False)

# End timer
end_time = time.time()
elapsed_time = end_time - start_time

# Convert elapsed time to hours and minutes
hours, rem = divmod(elapsed_time, 3600)
minutes, seconds = divmod(rem, 60)

print(f"RPG gibbs correction completed in {int(hours):02} hours, {int(minutes):02} minutes, {int(seconds):02} seconds.")
print(separator + "\n")

def group_alignment(group_list):
from mrtrix3 import app, run, path, image, fsl, MRtrixError
import numpy as np

fsl_suffix = fsl.suffix()

terminal_width = shutil.get_terminal_size().columns
separator = "=" * terminal_width

print("\n" + separator)
print('... Rigidly aligning groups...')
start_time = time.time()

# 1) find the longest series and extract the b0
group_len = []
Expand Down Expand Up @@ -218,6 +276,17 @@ def group_alignment(group_list):
show=False)
run.command('mrconvert -force dwi_align_series_aligned.mif working.mif', show=False)

# End timer
end_time = time.time()
elapsed_time = end_time - start_time

# Convert elapsed time to hours and minutes
hours, rem = divmod(elapsed_time, 3600)
minutes, seconds = divmod(rem, 60)

print(f"rigid alignment completed in {int(hours):02} hours, {int(minutes):02} minutes, {int(seconds):02} seconds.")
print(separator + "\n")


def convert_ants_xform(mat, i):
import scipy.io as sio
Expand All @@ -242,6 +311,13 @@ def pre_eddy_ants_moco(dwi_metadata):
import ants
from mrtrix3 import app, run, path, image, fsl, MRtrixError
import numpy as np

terminal_width = shutil.get_terminal_size().columns
separator = "=" * terminal_width

print("\n" + separator)
print('... ANTS motion correction...')
start_time = time.time()

run.command('mrconvert -force -export_grad_fsl working.bvec working.bval working.mif working.nii', show=False)
nii = ants.image_read('working.nii')
Expand Down Expand Up @@ -280,6 +356,17 @@ def pre_eddy_ants_moco(dwi_metadata):
np.savetxt('working_antsmoco.bvec', dirs_rot.T, delimiter=' ', fmt='%4.10f')
run.command('mrconvert -force -fslgrad working_antsmoco.bvec working.bval working_antsmoco.nii working.mif', show=False)

# End timer
end_time = time.time()
elapsed_time = end_time - start_time

# Convert elapsed time to hours and minutes
hours, rem = divmod(elapsed_time, 3600)
minutes, seconds = divmod(rem, 60)

print(f"motion correction completed in {int(hours):02} hours, {int(minutes):02} minutes, {int(seconds):02} seconds.")
print(separator + "\n")

def run_pre_align(dwi_metadata):
from mrtrix3 import app, run, path, image, fsl, MRtrixError

Expand Down Expand Up @@ -308,7 +395,13 @@ def run_eddy(shell_table, dwi_metadata):
import os
import glob

print("...Eddy current, EPI, motion correction...")
terminal_width = shutil.get_terminal_size().columns
separator = "=" * terminal_width

print("\n" + separator)
print('... Eddy current, EPI, motion correction...')
start_time = time.time()

run.command('mrconvert -force -export_grad_fsl working.bvec working.bval working.mif working.nii', show=False)

eddyopts = '" --cnr_maps --repol --data_is_shelled "'
Expand Down Expand Up @@ -733,15 +826,31 @@ def run_eddy(shell_table, dwi_metadata):
raise MRtrixError("the eddy option must run alongside -rpe_header, -rpe_all, or -rpe_pair option")

run.command('mrconvert -force -fslgrad working.bvec working.bval dwiec.mif working.mif', show=False)
# End timer
end_time = time.time()
elapsed_time = end_time - start_time

# Convert elapsed time to hours and minutes
hours, rem = divmod(elapsed_time, 3600)
minutes, seconds = divmod(rem, 60)

print(f"Eddy current, EPI, motion correction completed in {int(hours):02} hours, {int(minutes):02} minutes, {int(seconds):02} seconds.")
print(separator + "\n")

def run_b1correct(dwi_metadata):
from mrtrix3 import run

DWInlist = dwi_metadata['dwi_list']
idxlist = dwi_metadata['idxlist']

terminal_width = shutil.get_terminal_size().columns
separator = "=" * terminal_width

print("\n" + separator)
print('... B1 correction...')
start_time = time.time()

# b1 bias field correction
print("...B1 correction...")
if len(DWInlist) == 1:
run.command('dwibiascorrect ants -bias biasfield.mif working.mif dwibc.mif')
else:
Expand All @@ -759,6 +868,17 @@ def run_b1correct(dwi_metadata):
run.command('mrcat -axis 3 ' + DWImif + ' dwibc.mif')
run.command('mrconvert -force dwibc.mif working.mif', show=False)

# End timer
end_time = time.time()
elapsed_time = end_time - start_time

# Convert elapsed time to hours and minutes
hours, rem = divmod(elapsed_time, 3600)
minutes, seconds = divmod(rem, 60)

print(f"Eddy current, EPI, motion correction completed in {int(hours):02} hours, {int(minutes):02} minutes, {int(seconds):02} seconds.")
print(separator + "\n")

def create_brainmask(fsl_suffix):
import os, gzip, shutil
from mrtrix3 import run
Expand Down
33 changes: 31 additions & 2 deletions lib/designer_input_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -68,6 +68,21 @@ def get_input_info(input, fslbval, fslbvec, bids):
UserBidspath = bids.resplit(',')
bidslist = [os.path.realpath(i) for i in UserBidspath]

bshapelist = [i + '.bshape' for i in DWInlist]
telist = [i + '.echotime' for i in DWInlist]
try:
for i in bshapelist:
if not os.path.exists(i):
bshapelist = None
break
for i in telist:
if not os.path.exists(i):
telist = None
break
except:
bshapelist = None
telist = None

# if the user provides the path to phases
try:
phase_cpath = app.ARGS.phase.rsplit(',')
Expand All @@ -87,7 +102,10 @@ def get_input_info(input, fslbval, fslbvec, bids):
'dwi_list': DWInlist,
'phase_list': phase_nlist,
'phase_ext': phase_ext,
'bidslist': bidslist}
'bidslist': bidslist,
'bshapelist': bshapelist,
'telist': telist
}

return dwi_metadata

Expand Down Expand Up @@ -205,6 +223,10 @@ def assert_inputs(dwi_metadata, args_pe_dir, args_pf):
else:
TE = TE_app

#if TE is not None:
if (len(set(TE)) > 1) and (not app.ARGS.rpe_te):
raise MRtrixError('If data has variable echo time and no RPE TE is specified, please use the -rpe_te flag to specify the RPE TE')

# if no partial fourier information is found, assume full sampling
if not args_pf and not pf_bids:
args_pf = 1
Expand Down Expand Up @@ -268,6 +290,8 @@ def convert_input_data(dwi_metadata):
bshape_per_series = dwi_metadata['bshape']
phase_n_list = dwi_metadata['phase_list']
phase_ext = dwi_metadata['phase_ext']
bshapelist_input = dwi_metadata['bshapelist']
telist_input = dwi_metadata['telist']

if len(dwi_n_list) == 1:
if not isdicom:
Expand Down Expand Up @@ -312,7 +336,7 @@ def convert_input_data(dwi_metadata):
DWImif = ' '.join(miflist)
cmd = ('mrcat -axis 3 %s %s/dwi.mif' % (DWImif, app.SCRATCH_DIR))
run.command(cmd)

try:
if phase_n_list:
if len(phase_n_list) == 1:
Expand Down Expand Up @@ -356,6 +380,11 @@ def convert_input_data(dwi_metadata):
telist = [item for sublist in telist for item in sublist]
bshapelist = [item for sublist in bshapelist for item in sublist]

if telist_input is not None:
telist = np.hstack([np.loadtxt(i) for i in telist_input])
if bshapelist_input is not None:
bshapelist = np.hstack([np.loadtxt(i) for i in bshapelist_input])

dwi_metadata['idxlist'] = idxlist
dwi_metadata['echo_time_per_volume'] = np.array(telist)
dwi_metadata['bshape_per_volume'] = np.array(bshapelist)
Expand Down
Loading

0 comments on commit 2cc4ae2

Please sign in to comment.