invert_elastix

#!/bin/bash

SCRIPT_DIR="$( cd "$( dirname "${BASH_SOURCE[0]}" )" && pwd )"

# For compatibility with systems where realpath is not installed, e.g. MacOS.
if ! which realpath &> /dev/null; then
    realpath() {
        path="$1"
        # First, if the file is a symlink, recursively resolve it
        while [ -L "$path" ]; do
            link_path=$(readlink $path)
            case ${link_path:0:1} in
                "/") path="$link_path" ;;
                "~") path="$link_path" ;;
                *) path="$(dirname "$path")/$link_path"
            esac
        done
        # Then resolve folder symlinks and ".." and other junk in the full path
        echo "$(cd "$(dirname "$path")" && pwd -P)/$(basename "$1")"
    }
fi

# Load run_elastix settings
source run_elastix_settings

# Configure elastix path
if which elastix > /dev/null; then
    elastix_path=$(which elastix)
    echo "Using elastix found on the PATH: $elastix_path"
else
    elastix_path=$(ls $elastix_installation_location/elastix*/bin/elastix)
    n_elastix_found=$(ls $elastix_installation_location/elastix*/bin/elastix | wc -l)
    if [ "$n_elastix_found" -eq 0 ]; then
        >&2 echo "elastix could not be found at $elastix_installation_location/elastix*/bin/elastix"
        >&2 echo "Make sure you have downloaded elastix and then edited ${SCRIPT_DIR}/run_elastix_settings to specify where you put the elastix files"
        >&2 echo "See README.md for more details"
        exit 1
    elif [ "$n_elastix_found" -gt 1 ]; then
        >&2 echo "$n_elastix_found elastix executables found in $elastix_installation_location:"
        >&2 echo $elastix_path
        elastix_path=$(ls -t $elastix_installation_location/elastix*/bin/elastix | head -n 1)
        >&2 echo "Using the version with the most recent timestamp: $elastix_path"
    else
        >&2 echo "Using elastix found at: $elastix_path"
    fi
    # Add the folder containing elastix to the path
    export PATH=$(dirname "$elastix_path"):$PATH
fi

# Due to MacOS System Integrity Protection purging LD_* environment variables,
# we set it explicitly here so that elastix can find is library file
# libANNlib-5.0. Doing this also makes sure elastix will work on other OSes
# even if the user hasn't set their LD_LIBRARY_PATH to point to the folder
# containing elastix's libANNlib-5.0 
elastix_path=$(realpath $elastix_path)
# Allow elastix to search for libANNlib-5.0 in the elastix folder and in ../lib
export LD_LIBRARY_PATH=$(dirname "$elastix_path"):$(dirname "$(dirname "$elastix_path")")/lib


show_help () {
    >&2 echo "Usage: invert_elastix folder_containing_transform_to_invert [-v|--verbose] [-f|--fake]"
    >&2 echo "                      [-o output_suffix] [-n n_threads] [-s spacings] [-w|--write_intermediate_output]"
    >&2 echo "Example: ./invert_elastix 21G01LexA_vnc2_60x_T1_2019_09_19_neuropil_elastix_to_moving_template/elastix_Bspline/4spacing_20bendingweight"
}

if [ "$#" -eq 0 ] || [ "$1" == "--help" ] || [ "$1" == "-h" ]; then
    show_help
    exit 1
fi

if [ ! -f ${SCRIPT_DIR}/elastixParams_invertBspline.txt ]; then
    show_help
    >&2 echo "ERROR: elastixParams_invertBspline.txt must exist in ${SCRIPT_DIR}"
    exit 1
fi

# Check whether the GNU sed or BSD/MacOS sed is installed and use
# the form of the -i argument required by the detected version
# Information from https://stackoverflow.com/a/44864004
# Implementation from https://stackoverflow.com/a/9323261
if sed --version &> /dev/null; then
    dashi=(-i)  # GNU sed
else
    dashi=(-i "")  # BSD/MacOS sed
fi


#DEFAULTS FOR ARGUMENT-MODIFIABLE VARIABLES GO HERE
verbose=false
fake=false
output_suffix=""
moving_template=false
run_affine=true
run_Bspline=true
n_threads=$number_of_threads_to_use
spacings="6"
write_intermediate_output=false


positionalArgs=()
unknownOptions=()
while [ "$#" -gt 0 ]; do
    case "$1" in
        -v)
            verbose=true
            >&2 echo "Running verbose"
            shift
        ;;
        -f)
            fake=true
            >&2 echo "fake=true"
            shift
        ;;
        -o)
            output_suffix="$2"
            shift; shift
        ;;
        -n)
            n_threads="$2"
            shift; shift
        ;;
        -s)
            spacings="$2"
            shift; shift
        ;;
        -w)
            write_intermediate_output=true
            shift
        ;;
        #-e|--exampleLongOption)
        #    thisOption=true
        #    shift
        #;;
        #-n|--number) #For flags that take a separate argument (like "-n 8" or any other "-[letter] [argument]"), use this template
        #    number="$2"
        #    shift; shift
        #;;
        *)
            if [ " ${1:0:1}" == " -" ]; then #Ignore arguments starting with - that aren't explicitly listed above
                if [ "$#" -eq 1 ]; then
                    show_help
                    >&2 echo "ERROR: The only parameter is an unknown option. Aborting."
                    exit 1
                fi
                unknownOptions+=("$1")
                >&2 echo "WARNING: Unknown option $1, ignoring"
            else
                positionalArgs+=("$1") #Store arguments (other than ones recognized above) in order
            fi
            shift
        ;;
    esac
done
set -- "${positionalArgs[@]}"
if [ "$#" -ne 1 ]; then
    >&2 echo "Provide exactly 1 positional argument"
    show_help
    exit 1
fi
#echo $@

#FAKEABLE COMMANDS GO HERE.
fakeableCommands="rm mv cp ln chmod mkdir sed exit elastix transformix"
#For the commands listed here, you can use "$command" instead of "command" anywhere in your code that you want the command to echo instead of execute when running in fake mode. For instance, write the line "$mv file1 file2" instead of "mv file1 file2" if you want that line's mv command to not execute when the -f flag is given. If you write "mv file1 file2", then your mv command will happen normally regardless of if -f is given.

if $fake; then
    for cmd in $fakeableCommands; do
        eval $cmd=\"echo \$cmd\"
    done
else
    for cmd in $fakeableCommands; do
        eval $cmd=\"\$cmd\"
    done
fi


#CODE GOES HERE

base=$1 #$(realpath "$1")
if [ ! -d "$base" ]; then
    show_help
    >&2 echo "ERROR: Argument must be an existing folder"
    exit 1
fi

elastix_log=$base/elastix.log
if [ ! -f "$elastix_log" ]; then
    show_help
    >&2 echo "ERROR: elastix.log not present in $base"
    exit 1
fi

moving=$(grep -- "-m  " "$elastix_log" | sed 's/-m *//')
if [ ! -f "$moving" ]; then
    moving=$(realpath $base/../../../$moving)
    if [ ! -f  "$moving" ]; then
        show_help
        >&2 echo "ERROR: moving image could not be found at $moving"
        exit 1
    fi
fi
fixed=$(grep -- "-f  " "$elastix_log" | sed 's/-f *//')
if [ ! -f "$fixed" ]; then
    fixed=$(realpath "$base/../../../$fixed")
    if [ ! -f  "$fixed" ]; then
        show_help
        >&2 echo "ERROR: fixed image could not be found at $fixed"
        exit 1
    fi
fi

echo "Fixed: $fixed"
echo "Moving: $moving"

if [ ! -f "$base/TransformParameters.0.txt" ]; then
    show_help
    echo "ERROR: TransformationParameters.0.txt not present in $base"
    exit 1
fi

if [ -n "$output_suffix" ]; then
    output_suffix=_$output_suffix #Put a _ in just before the suffix if it's non-empty
fi

for spacing in $spacings; do
    output_folder=$base/inverted_${spacing}spacing$output_suffix

    mkdir -p $output_folder
    elastixParams="$output_folder/elastixParams_invertBspline.txt"

    if [ -e "$elastixParams" -a ! -w "$elastixParams" ]; then
        >&2 read -p "Parameter file already present ($elastixParams) but not overwritable. Are you sure you want to overwrite it? (y/n) " input
        case $input in [yY]) echo "Overwriting"; chmod +w "$elastixParams";; *) echo "Exiting"; exit 1 ;; esac
    fi 
    cp "${SCRIPT_DIR}/elastixParams_invertBspline.txt" "$elastixParams"
    sed "${dashi[@]}" "s/(FinalGridSpacingInPhysicalUnits _GRIDSPACING_)/(FinalGridSpacingInPhysicalUnits ${spacing//-/ })/" "$elastixParams"
    if $write_intermediate_output; then sed "${dashi[@]}" 's/(WriteResultImage "false")/(WriteResultImage "true")/' "$elastixParams"; fi
    $chmod -w "$elastixParams"

    echo "Inverting $base/TransformParameters.0.txt with GridSpacing $spacing"
    #To invert a transform, -f and -m are both set to be the original moving image, and -t0 is set to the transform.
    #This is opposite of what the elastix manual says to do (where it says to set -f and -m to the fixed image)
    #but I tested both ways, and the inversion only worked correctly for me with -f and -m set to the moving image.
    $elastix -f "$moving" -m "$moving" -out "$output_folder/" \
        -t0 "$base/TransformParameters.0.txt" -threads $n_threads -p "$elastixParams"

    if [ ! -f "$output_folder/TransformParameters.0.txt" ]; then
        show_help
        echo "ERROR: Elastix failed, stopping invert_elastix"
        $exit 1
    fi

    if [ -f "$output_folder/result.0.nrrd" ]; then $mv "$output_folder/result.0.nrrd" "$output_folder/result.0.composedForwardAndInverseTransforms.nrrd"; fi

    $cp "$output_folder/TransformParameters.0.txt" "$output_folder/TransformParameters.1.txt"
    $sed "${dashi[@]}" 's/(InitialTransformParametersFileName.*)/(InitialTransformParametersFileName "NoInitialTransform")/' "$output_folder/TransformParameters.1.txt"

    #Apply the inverted transform to the fixed image
    $transformix -in "$fixed" \
        -out "$output_folder/" \
        -tp "$output_folder/TransformParameters.1.txt" \
        -threads $n_threads

    if [ ! -f "$output_folder/result.nrrd" ]; then
        show_help
        >&2 echo "ERROR: Transformix failed. Exiting."
        exit 1
    fi

    $mv "$output_folder/result.nrrd" "$output_folder/result.1.nrrd"
done