python_c_extension/_disptools.c
/*!
 * Method to compute a displacement with given Jacobian.
 */
Method to compute a displacement with given Jacobian.
Method to compute a displacement with given Jacobian .
static PyObject *method_displacement(PyObject *self, PyObject *args)
{
    FLOATING spacing[3];
    PyObject *spacing_tuple = NULL;
    PyArrayObject *field = NULL, *jacobian = NULL, *mask = NULL;
    FLOATING epsilon, tolerance, eta, eta_max, alpha, beta, gamma, delta, zeta, theta, iota;
    bool strict;
    long it_max;
    const char *algorithm;
    int gpu_id;
    (void) self;
    if (!PyArg_ParseTuple(args,
                          PY_FMT_DISPLACEMENT,
                          &spacing_tuple,
                          &jacobian,
                          &mask,
                          &epsilon,
                          &tolerance,
                          &eta,
                          &eta_max,
                          &alpha,
                          &beta,
                          &gamma,
                          &delta,
                          &zeta,
                          &theta,
                          &iota,
                          &strict,
                          &it_max,
                          &field,
                          &algorithm,
                          &gpu_id)) {
        return NULL;
    }
    if (strcmp(algorithm, ALGORITHM_GRADIENT) &&
        strcmp(algorithm, ALGORITHM_GREEDY) &&
        strcmp(algorithm, ALGORITHM_MATCHING))
    {
        PyErr_Format(PyExc_ValueError, "Invalid algorithm '%s'.", algorithm);
        return NULL;
    }
    if (!jacobian || !mask || !field || !spacing_tuple) {
        return NULL;
    }
    if (read_triplet(spacing_tuple, spacing, 3, "spacing")) {
        return NULL;
    }
    const int nz = PyArray_DIM(jacobian, 0);
    const int ny = PyArray_DIM(jacobian, 1);
    const int nx = PyArray_DIM(jacobian, 2);
    FLOATING *jacobian_data = (FLOATING*) PyArray_DATA(jacobian);
    bool *mask_data = (bool*) PyArray_DATA(mask);
    FLOATING *field_data = (FLOATING*) PyArray_DATA(field);
    DisplacementFunction displacement_function = NULL;
    if (!DISPTOOLS_HAS_CUDA && gpu_id >= 0) {
        PyErr_Format(PyExc_ValueError,
                     "Parameter gpu_id is set to %d, but this version "
                     "of disptools is built with no CUDA support.",
                     gpu_id);
        return NULL;
    }
    else if (DISPTOOLS_HAS_CUDA && gpu_id >= 0) {
#if DISPTOOLS_HAS_CUDA
        if(!set_device(gpu_id)) {
            PyErr_Format(PyExc_ValueError, "Cannot set CUDA device %d", gpu_id);
            return NULL;
        }
        if (!strcmp(algorithm, ALGORITHM_GRADIENT)) {
            displacement_function = &generate_displacement_gradient_cuda;
        }
        else if (!strcmp(algorithm, ALGORITHM_GREEDY)) {
            displacement_function = &generate_displacement_greedy_cuda;
        }
        else if (!strcmp(algorithm, ALGORITHM_MATCHING)) {
            PyErr_Format(PyExc_ValueError,
                         "Algorithm '%s' not supported on GPU",
                         algorithm);
            return NULL;
        }
#endif 
    }
    else {
        if (!strcmp(algorithm, ALGORITHM_GRADIENT)) {
            displacement_function = &generate_displacement_gradient;
        }
        else if (!strcmp(algorithm, ALGORITHM_GREEDY)) {
            displacement_function = &generate_displacement_greedy;
        }
        else if (!strcmp(algorithm, ALGORITHM_MATCHING)) {
            displacement_function = &volume_matching_3d;
        }
    }
    displacement_function(
            nx, ny, nz,
            spacing[0], spacing[1], spacing[2],
            (void*) jacobian_data,
            (void*) mask_data,
            epsilon,
            tolerance,
            eta,
            eta_max,
            alpha,
            beta,
            gamma,
            delta,
            zeta,
            theta,
            iota,
            strict,
            it_max,
            (void*) field_data);
    if (disptools_has_error()) {
        disptools_error_state *err = disptools_get_error();
        PyErr_Format(PyExc_RuntimeError,
                     "%s (%s:%d): %s\n%s",
                     err->function,
                     err->file,
                     err->line,
                     err->message,
                     err->trace);
        return NULL;
    }
    return Py_None;
}
static PyObject * method_displacement ( PyObject * self , PyObject * args ) { FLOATING spacing [ 3 ] ; PyObject * spacing_tuple = NULL ; PyArrayObject * field = NULL , * jacobian = NULL , * mask = NULL ; FLOATING epsilon , tolerance , eta , eta_max , alpha , beta , gamma , delta , zeta , theta , iota ; bool strict ; long it_max ; const char * algorithm ; int gpu_id ; ( void ) self ; if ( ! PyArg_ParseTuple ( args , PY_FMT_DISPLACEMENT , & spacing_tuple , & jacobian , & mask , & epsilon , & tolerance , & eta , & eta_max , & alpha , & beta , & gamma , & delta , & zeta , & theta , & iota , & strict , & it_max , & field , & algorithm , & gpu_id ) ) { return NULL ; } if ( strcmp ( algorithm , ALGORITHM_GRADIENT ) && strcmp ( algorithm , ALGORITHM_GREEDY ) && strcmp ( algorithm , ALGORITHM_MATCHING ) ) { PyErr_Format ( PyExc_ValueError , " " , algorithm ) ; return NULL ; } if ( ! jacobian || ! mask || ! field || ! spacing_tuple ) { return NULL ; } if ( read_triplet ( spacing_tuple , spacing , 3 , " " ) ) { return NULL ; } const int nz = PyArray_DIM ( jacobian , 0 ) ; const int ny = PyArray_DIM ( jacobian , 1 ) ; const int nx = PyArray_DIM ( jacobian , 2 ) ; FLOATING * jacobian_data = ( FLOATING * ) PyArray_DATA ( jacobian ) ; bool * mask_data = ( bool * ) PyArray_DATA ( mask ) ; FLOATING * field_data = ( FLOATING * ) PyArray_DATA ( field ) ; DisplacementFunction displacement_function = NULL ; if ( ! DISPTOOLS_HAS_CUDA && gpu_id >= 0 ) { PyErr_Format ( PyExc_ValueError , " " " " , gpu_id ) ; return NULL ; } else if ( DISPTOOLS_HAS_CUDA && gpu_id >= 0 ) { #if DISPTOOLS_HAS_CUDA if ( ! set_device ( gpu_id ) ) { PyErr_Format ( PyExc_ValueError , " " , gpu_id ) ; return NULL ; } if ( ! strcmp ( algorithm , ALGORITHM_GRADIENT ) ) { displacement_function = & generate_displacement_gradient_cuda ; } else if ( ! strcmp ( algorithm , ALGORITHM_GREEDY ) ) { displacement_function = & generate_displacement_greedy_cuda ; } else if ( ! strcmp ( algorithm , ALGORITHM_MATCHING ) ) { PyErr_Format ( PyExc_ValueError , " " , algorithm ) ; return NULL ; } #endif } else { if ( ! strcmp ( algorithm , ALGORITHM_GRADIENT ) ) { displacement_function = & generate_displacement_gradient ; } else if ( ! strcmp ( algorithm , ALGORITHM_GREEDY ) ) { displacement_function = & generate_displacement_greedy ; } else if ( ! strcmp ( algorithm , ALGORITHM_MATCHING ) ) { displacement_function = & volume_matching_3d ; } } displacement_function ( nx , ny , nz , spacing [ 0 ] , spacing [ 1 ] , spacing [ 2 ] , ( void * ) jacobian_data , ( void * ) mask_data , epsilon , tolerance , eta , eta_max , alpha , beta , gamma , delta , zeta , theta , iota , strict , it_max , ( void * ) field_data ) ; if ( disptools_has_error ( ) ) { disptools_error_state * err = disptools_get_error ( ) ; PyErr_Format ( PyExc_RuntimeError , " \n " , err -> function , err -> file , err -> line , err -> message , err -> trace ) ; return NULL ; } return Py_None ; }
Method to compute a displacement with given Jacobian.
Method to compute a displacement with given Jacobian .
/* Parse arguments */ /* Sanity check */ /* Error check */ /* Get the spacing */ /* Get the shape */ /* Get the data pointers */ /* Select the appropriate function for the algorithm */ // DISPTOOLS_HAS_CUDA /* Call the library function to compute the displacement */