diff --git a/README.md b/README.md index e5d5b8d..018def7 100644 --- a/README.md +++ b/README.md @@ -18,6 +18,7 @@ Takes TIF images (tiled or prestitched) and outputs a validated BIDS Microscopy - Python >= 3.11 - Lightsheet data: - Raw Ultramicroscope Blaze OME TIFF files (include `blaze` in the acquisition tag) + - can be 2D or 3D TIFF files - Prestitched TIFF files (include `prestitched` in the acquisition tag) @@ -63,10 +64,8 @@ or for snakemake<8.0, use: snakemake -c all --use-singularity ``` -Note: if you run the workflow on a system with large memory, you will need to set the heap size for the stitching and fusion rules. This can be done with e.g.: `--set-resources bigstitcher_spark_stitching:mem_mb=60000 bigstitcher_spark_fusion:mem_mb=100000` +Note: if you run the workflow on a system with large memory, you will need to set the heap size for the stitching and fusion rules. This can be done with e.g.: `--set-resources bigstitcher_stitching:mem_mb=60000 bigstitcher_fusion:mem_mb=100000` 7. If you want to run the workflow using a batch job submission server, please see the executor plugins here: https://snakemake.github.io/snakemake-plugin-catalog/ - -Alternate usage of this workflow (making use of conda) is described in the [Snakemake Workflow Catalog](https://snakemake.github.io/snakemake-workflow-catalog?repo=khanlab/SPIMprep). diff --git a/config/config.yml b/config/config.yml index edc16ba..ae096d2 100644 --- a/config/config.yml +++ b/config/config.yml @@ -36,12 +36,10 @@ bigstitcher: downsample_in_x: 4 downsample_in_y: 4 downsample_in_z: 1 - method: "phase_corr" #unused - methods: #unused + methods: #unused, only for reference phase_corr: "Phase Correlation" optical_flow: "Lucas-Kanade" filter_pairwise_shifts: - enabled: 1 #unused min_r: 0.7 max_shift_total: 50 global_optimization: @@ -64,7 +62,7 @@ bigstitcher: block_size_factor_z: 32 ome_zarr: - desc: sparkstitchedflatcorr + desc: stitchedflatcorr max_downsampling_layers: 5 # e.g. 4 levels: { 0: orig, 1: ds2, 2: ds4, 3: ds8, 4: ds16} rechunk_size: #z, y, x - 1 @@ -155,5 +153,5 @@ report: containers: - spimprep: 'docker://khanlab/spimprep-deps:main' + spimprep: 'docker://khanlab/spimprep-deps:v0.1.0' diff --git a/workflow/macros/AutostitchMacro.ijm b/workflow/macros/AutostitchMacro.ijm deleted file mode 100644 index e9df4df..0000000 --- a/workflow/macros/AutostitchMacro.ijm +++ /dev/null @@ -1,56 +0,0 @@ -args = getArgument() -args = split(args, " "); - -dataset_xml = args[0] -pairwise_method = args[1] -ds_x = args[2] -ds_y = args[3] -ds_z = args[4] -do_filter = args[5] -min_r = args[6] -do_global = args[7] -global_strategy = args[8] - -run("Calculate pairwise shifts ...", - "select=" + dataset_xml + - " process_angle=[All angles] " + - " process_channel=[All channels] " + - " process_illumination=[All illuminations] " + - " process_tile=[All tiles] " + - " process_timepoint=[All Timepoints] " + - " method=["+ pairwise_method +"] " + - " channels=[Average Channels] " + - " downsample_in_x=" + ds_x + - " downsample_in_y=" + ds_y + - " downsample_in_z=" + ds_z ); - - -if ( do_filter == 1 ){ - -run("Filter pairwise shifts ...", - "select=" + dataset_xml + - " min_r=" + min_r + - " max_r=1 " + - " max_shift_in_x=0 " + - " max_shift_in_y=0 " + - " max_shift_in_z=0 " + - " max_displacement=0"); -} - -if ( do_global == 1 ){ - -run("Optimize globally and apply shifts ...", - "select=" + dataset_xml + - " process_angle=[All angles] " + - " process_channel=[All channels] " + - " process_illumination=[All illuminations] " + - " process_tile=[All tiles] " + - " process_timepoint=[All Timepoints] " + - " relative=2.500 " + - " absolute=3.500 " + - " global_optimization_strategy=["+global_strategy+"] " + - " fix_group_0-0,"); -} - -// quit after we are finished -eval("script", "System.exit(0);"); diff --git a/workflow/macros/FuseImageMacroZarr.ijm b/workflow/macros/FuseImageMacroZarr.ijm deleted file mode 100644 index 9d49b06..0000000 --- a/workflow/macros/FuseImageMacroZarr.ijm +++ /dev/null @@ -1,47 +0,0 @@ -args = getArgument() -args = split(args, " "); - -in_dataset_xml = args[0] -downsampling = args[1] -channel = args[2] -output_zarr = args[3] -block_size_x = args[4] -block_size_y = args[5] -block_size_z = args[6] -block_size_factor_x = args[7] -block_size_factor_y = args[8] -block_size_factor_z = args[9] - -run("Fuse dataset ...", - "select=" + in_dataset_xml + - " process_angle=[All angles] " + - " process_channel=[Single channel (Select from List)] " + - " processing_channel=[channel " + channel + "]" + - " process_illumination=[All illuminations] " + - " process_tile=[All tiles] " + - " process_timepoint=[All Timepoints] " + - " bounding_box=[Currently Selected Views] " + - " preserve_original " + - " downsampling=" + downsampling + - " interpolation=[Linear Interpolation] " + - " pixel_type=[16-bit unsigned integer] " + - " interest_points_for_non_rigid=[-= Disable Non-Rigid =-] " + - " blend produce=[Each timepoint & channel] " + - " fused_image=[ZARR/N5/HDF5 export using N5-API] " + - " define_input=[Auto-load from input data (values shown below)] " + - " export=ZARR create_0 " + - " zarr_dataset_path=" + output_zarr + - " zarr_base_dataset=/ " + - " zarr_dataset_extension=/s0 " + - " show_advanced_block_size_options " + - " block_size_x=" + block_size_x + - " block_size_y=" + block_size_y + - " block_size_z=" + block_size_z + - " block_size_factor_x=" + block_size_factor_x + - " block_size_factor_y=" + block_size_factor_y + - " block_size_factor_z=" + block_size_factor_z + - " subsampling_factors=[{ {1,1,1}}]"); - - -// quit after we are finished -eval("script", "System.exit(0);"); diff --git a/workflow/rules/bigstitcher.smk b/workflow/rules/bigstitcher.smk index 4256e3b..24a4a3e 100644 --- a/workflow/rules/bigstitcher.smk +++ b/workflow/rules/bigstitcher.smk @@ -94,74 +94,7 @@ rule zarr_to_bdv: "../scripts/zarr_to_n5_bdv.py" -rule bigstitcher: - input: - dataset_n5=rules.zarr_to_bdv.output.bdv_n5, - dataset_xml=rules.zarr_to_bdv.output.bdv_xml, - ijm=Path(workflow.basedir) / "macros" / "AutostitchMacro.ijm", - params: - fiji_launcher_cmd=get_fiji_launcher_cmd, - macro_args=get_macro_args_bigstitcher, - rm_old_xml=lambda wildcards, output: f"rm -f {output.dataset_xml}~?", - output: - launcher=temp( - bids( - root=work, - subject="{subject}", - datatype="micr", - sample="{sample}", - acq="{acq}", - desc="{desc}", - suffix="bigstitcherfiji.sh", - ) - ), - dataset_xml=temp( - bids( - root=work, - subject="{subject}", - datatype="micr", - sample="{sample}", - acq="{acq}", - desc="{desc}", - suffix="bigstitcherfiji.xml", - ) - ), - benchmark: - bids( - root="benchmarks", - datatype="bigstitcherfiji", - subject="{subject}", - sample="{sample}", - acq="{acq}", - desc="{desc}", - suffix="benchmark.tsv", - ) - log: - bids( - root="logs", - datatype="bigstitcherfiji", - subject="{subject}", - sample="{sample}", - acq="{acq}", - desc="{desc}", - suffix="log.txt", - ), - container: - config["containers"]["spimprep"] - resources: - runtime=30, - mem_mb=10000, - threads: config["cores_per_rule"] - group: - "preproc" - shell: - "cp {input.dataset_xml} {output.dataset_xml} && " - " {params.fiji_launcher_cmd} && " - " echo ' -macro {input.ijm} \"{params.macro_args}\"' >> {output.launcher} " - " && {output.launcher} |& tee {log} && {params.rm_old_xml}" - - -rule bigstitcher_spark_stitching: +rule bigstitcher_stitching: input: dataset_n5=rules.zarr_to_bdv.output.bdv_n5, dataset_xml=rules.zarr_to_bdv.output.bdv_xml, @@ -230,10 +163,10 @@ rule bigstitcher_spark_stitching: "{params.rm_old_xml}" -rule bigstitcher_spark_solver: +rule bigstitcher_solver: input: dataset_n5=rules.zarr_to_bdv.output.bdv_n5, - dataset_xml=rules.bigstitcher_spark_stitching.output.dataset_xml, + dataset_xml=rules.bigstitcher_stitching.output.dataset_xml, params: downsampling="--downsampling={dsx},{dsy},{dsz}".format( dsx=config["bigstitcher"]["calc_pairwise_shifts"]["downsample_in_x"], @@ -295,7 +228,6 @@ rule bigstitcher_spark_solver: "{params.rm_old_xml}" #lambda 0.1 is default (can expose this if needed) - rule bigstitcher_fusion: input: dataset_n5=bids( @@ -320,98 +252,6 @@ rule bigstitcher_fusion: else "stitching" ), ), - ijm=Path(workflow.basedir) / "macros" / "FuseImageMacroZarr.ijm", - params: - fiji_launcher_cmd=get_fiji_launcher_cmd, - macro_args=get_macro_args_zarr_fusion, - output: - launcher=temp( - bids( - root=work, - subject="{subject}", - datatype="micr", - sample="{sample}", - acq="{acq}", - desc="stitched{desc}", - stain="{stain}", - suffix="fuseimagen5.sh", - ) - ), - zarr=temp( - directory( - bids( - root=work, - subject="{subject}", - datatype="micr", - sample="{sample}", - acq="{acq}", - desc="stitched{desc}", - stain="{stain}", - suffix="SPIM.zarr", - ) - ) - ), - benchmark: - bids( - root="benchmarks", - datatype="fuse_dataset", - subject="{subject}", - sample="{sample}", - acq="{acq}", - desc="stitched{desc}", - stain="{stain}", - suffix="benchmark.tsv", - ) - log: - bids( - root="logs", - datatype="fuse_dataset", - subject="{subject}", - sample="{sample}", - acq="{acq}", - desc="stitched{desc}", - stain="{stain}", - suffix="log.txt", - ), - container: - config["containers"]["spimprep"] - resources: - runtime=30, - mem_mb=40000, - threads: config["cores_per_rule"] - group: - "preproc" - shell: - " {params.fiji_launcher_cmd} && " - " echo ' -macro {input.ijm} \"{params.macro_args}\"' >> {output.launcher} " - " && {output.launcher} |& tee {log}" - - -rule bigstitcher_spark_fusion: - input: - dataset_n5=bids( - root=work, - subject="{subject}", - datatype="micr", - sample="{sample}", - acq="{acq}", - desc="{desc}", - suffix="bdv.n5", - ), - dataset_xml=bids( - root=work, - subject="{subject}", - datatype="micr", - sample="{sample}", - acq="{acq}", - desc="{desc}", - suffix="bigstitcher{}.xml".format( - "solver" - if config["bigstitcher"]["global_optimization"]["enabled"] - else "stitching" - ), - ), - ijm=Path(workflow.basedir) / "macros" / "FuseImageMacroZarr.ijm", params: channel=lambda wildcards: "--channelId={channel}".format( channel=get_stains(wildcards).index(wildcards.stain) @@ -438,7 +278,7 @@ rule bigstitcher_spark_fusion: datatype="micr", sample="{sample}", acq="{acq}", - desc="sparkstitched{desc}", + desc="stitched{desc}", stain="{stain}", suffix="SPIM.zarr", ) diff --git a/workflow/rules/common.smk b/workflow/rules/common.smk index 88d40e7..80a773a 100644 --- a/workflow/rules/common.smk +++ b/workflow/rules/common.smk @@ -257,52 +257,6 @@ def get_stains(wildcards): return df.iloc[0][stain_columns].dropna().tolist() -# bigstitcher -def get_fiji_launcher_cmd(wildcards, output, threads, resources): - launcher_opts_find = "-Xincgc" - launcher_opts_replace = f"-XX:+UseG1GC -verbose:gc -XX:+PrintGCDateStamps -XX:ActiveProcessorCount={threads}" - pipe_cmds = [] - pipe_cmds.append("ImageJ-linux64 --dry-run --headless --console") - pipe_cmds.append(f"sed 's/{launcher_opts_find}/{launcher_opts_replace}'/") - pipe_cmds.append( - rf"sed 's/-Xmx[0-9a-z]\+/-Xmx{resources.mem_mb}m -Xms{resources.mem_mb}m/'" - ) - pipe_cmds.append("tr --delete '\\n'") - return "|".join(pipe_cmds) + f" > {output.launcher} && chmod a+x {output.launcher} " - - -def get_macro_args_bigstitcher(wildcards, input, output): - return "{dataset_xml} {pairwise_method} {ds_x} {ds_y} {ds_z} {do_filter} {min_r} {do_global} {global_strategy}".format( - dataset_xml=output.dataset_xml, - pairwise_method=config["bigstitcher"]["calc_pairwise_shifts"]["methods"][ - config["bigstitcher"]["calc_pairwise_shifts"]["method"] - ], - ds_x=config["bigstitcher"]["calc_pairwise_shifts"]["downsample_in_x"], - ds_y=config["bigstitcher"]["calc_pairwise_shifts"]["downsample_in_y"], - ds_z=config["bigstitcher"]["calc_pairwise_shifts"]["downsample_in_z"], - do_filter=config["bigstitcher"]["filter_pairwise_shifts"]["enabled"], - min_r=config["bigstitcher"]["filter_pairwise_shifts"]["min_r"], - do_global=config["bigstitcher"]["global_optimization"]["enabled"], - global_strategy=config["bigstitcher"]["global_optimization"]["strategies"][ - config["bigstitcher"]["global_optimization"]["strategy"] - ], - ) - - -def get_macro_args_zarr_fusion(wildcards, input, output): - return "{dataset_xml} {downsampling} {channel:02d} {output_zarr} {bsx} {bsy} {bsz} {bsfx} {bsfy} {bsfz}".format( - dataset_xml=input.dataset_xml, - downsampling=config["bigstitcher"]["fuse_dataset"]["downsampling"], - channel=get_stains(wildcards).index(wildcards.stain), - output_zarr=output.zarr, - bsx=config["bigstitcher"]["fuse_dataset"]["block_size_x"], - bsy=config["bigstitcher"]["fuse_dataset"]["block_size_y"], - bsz=config["bigstitcher"]["fuse_dataset"]["block_size_z"], - bsfx=config["bigstitcher"]["fuse_dataset"]["block_size_factor_x"], - bsfy=config["bigstitcher"]["fuse_dataset"]["block_size_factor_y"], - bsfz=config["bigstitcher"]["fuse_dataset"]["block_size_factor_z"], - ) - def get_output_ome_zarr_uri(): if is_remote(config["root"]):