Snakemake workflow: Reconstructing raw tomography data using tomopy.
Help improve this workflow!
This workflow has been published but could be further improved with some additional meta data:- Keyword(s) in categories input, output
You can help improve this workflow by suggesting the addition or removal of keywords, suggest changes and report issues, or request to become a maintainer of the Workflow .
Snakemake workflow: Reconstructing raw tomography data
A Snakemake worfklow for tomographically reconstructing raw data using tomopy .
Installation
First download this repo and navigate to it
git clone https://codebase.helmholtz.cloud/gernha62/reconstructing-raw-tomography-data.git
cd /path/to/repo
(Optional) Download the example folder with:
wget -m -np https://doi2.psi.ch/datasets/das/work/p15/p15869/compression/MI04_02/tif
Create a virtual environment and install all necessary packages (requires conda):
conda env create --name reconstr_env --file workflow/envs/reconstr.yml
Activate the new virtual environment:
conda activate reconstr_env
Configuration
To configure the workflow, adapt the config file found at
config/config.yaml
. The config looks as follows:
number_of_darks: 50
number_of_flats: 100
number_of_projections: 501
rotation_center: 508.77
raw_data: MI04_02: doi2.psi.ch/datasets/das/work/p15/p15869/compression/MI04_02/tif
In the config, adjust
number_of_darks
,
number_of_flats
,
number_of_projections
and
rotation_center
to the number of darks, flats, projections and the rotation center of your dataset. The necessary information can usually be found in the .log file of the folder that contains the raw data.
MI04_02: doi2.psi.ch/datasets/das/work/p15/p15869/compression/MI04_02/tif
denotes the path to the example folder used for reconstruction and the keyword
MI04_02
will be used to name the output (e.g. in this case the output folder will be named
recon_dir_MI04_02
). Replace the examle path with the path to the dataset you want to reconstruct. Additionally, if you want the name of the output folder to have a different suffix, replace the keyword
MI04_02
with a name you prefer.
Run the workflow
If the .tif files contain a numerical prefix that is not separated from the actual image index, it is best to first rename the files. The files will be renamed to
00001.tif
,
00002.tif
and so on. If the renaming is needed, run:
snakemake --cores 1 'logs/renamefile_MI04_02.log'
If you replaced the keyword
MI04_02
in the config file then adjust the command accordingly (e.g. if you replaced the keyword with
Tomo_dataset
then the command should be
snakemake --cores 1 'logs/renamefile_Tomo_dataset.log'
).
Before trying to compute the reconstructions, make sure you have enough memory available (ideally more than 60 GB). To compute the reconstructions using one core, use the command:
snakemake --cores 1
If you want to use all available cores instead, use:
snakemake --cores all
This creates a folder in
results
with the reconstructed data.
Code Snippets
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 | import tomopy import dxchange import numpy as np from glob import glob def reconstruction(path_to_output, tifFile, numDarks, numFlats, numProjs, rot_center_init): #Load Raw Data into Memory tifRef = tifFile firstIndex = 1 #Read Darks darkStart = firstIndex darkStop = darkStart + numDarks dark = dxchange.reader.read_tiff_stack(tifRef, [i for i in range(darkStart, darkStop)]) #Read Flats flatStart = darkStop flatStop = flatStart + numFlats flat = dxchange.reader.read_tiff_stack(tifRef, [i for i in range(flatStart, flatStop)]) #Read Projections projStart = flatStop projStop = projStart + numProjs proj = dxchange.reader.read_tiff_stack(tifRef, [i for i in range(projStart, projStop)]) print ("Shape of projection cube", proj.shape) #Define equiangular angles (angles here not part of raw data) theta = tomopy.angles(proj.shape[0]) #Start Reconstruction #Correct for dark and flat threads = 12 proj = tomopy.normalize(proj, flat, dark, ncore=threads) proj = tomopy.minus_log(proj, ncore=threads) rot_center = rot_center_init #Inverse Radon recon = tomopy.recon(proj, theta, center=rot_center, algorithm='gridrec', sinogram_order=False, ncore=threads) #Cosmetics recon = tomopy.circ_mask(recon, axis=0, ratio=0.95, ncore=threads) print("Reconstruction finished") dxchange.write_tiff_stack(recon, fname=path_to_output + '/recon') print ("Finished") def main(): path_to_directory = snakemake.input[0] fnames = glob(path_to_directory + '/' + '*.tif') tiff_file_example = fnames[0] path_to_output = snakemake.output[0] numDarks_magma = snakemake.params.numDarks numFlats_magma = snakemake.params.numFlats numProjs_magma = snakemake.params.numProjs rot_center_magma = snakemake.params.rotCenter reconstruction(path_to_output, tiff_file_example, numDarks_magma, numFlats_magma, numProjs_magma, rot_center_magma) if __name__ == '__main__': main() |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 | import os from glob import glob import sys #numerical prefix that is not separated from the actual image index confuses the dxchange filename parser, # so it is better to first rename the files. This function renames the files into 00000.tif, 00001.tif, ... def rename_tifs(path_to_directory): fnames = glob(path_to_directory + '/' + '*.tif') fnames.sort() for f in range(len(fnames)): new_name = path_to_directory + '/' + '{:0>5}'.format(str(f+1)) + '.tif' os.rename(fnames[f], new_name) with open(snakemake.log[0], "w") as f: sys.stderr = sys.stdout = f f.write('success') def main(): path_to_directory = snakemake.input[0] rename_tifs(path_to_directory) if __name__ == '__main__': main() |
20 21 | script: 'scripts/rename_tifs.py' |
33 34 | script: 'scripts/reconstructs_tomo_dataset.py' |
Support
- Future updates