Commit 4d9599e8 authored by Marc Modat's avatar Marc Modat

(1) When using OpenMP on windows, all the incremented variables are defined as...

(1) When using OpenMP on windows, all the incremented variables are defined as long. (2) reg_aladin now returns an error when the variance and inlier argument are outside of the range [0;100]
parent 1b3fc7e1
......@@ -30,7 +30,10 @@ set(NiftyReg_VERSION_MINOR 3)
set(NiftyReg_VERSION_PATCH 9)
#-----------------------------------------------------------------------------
if(COMMAND cmake_policy)
cmake_policy(SET CMP0003 NEW)
cmake_policy(SET CMP0003 NEW)
if("${CMAKE_MAJOR_VERSION}.${CMAKE_MINOR_VERSION}.${CMAKE_PATCH_VERSION}" MATCHES "^2\\.8\\.12$")
cmake_policy(SET CMP0022 OLD)
endif("${CMAKE_MAJOR_VERSION}.${CMAKE_MINOR_VERSION}.${CMAKE_PATCH_VERSION}" MATCHES "^2\\.8\\.12$")
endif(COMMAND cmake_policy)
#-----------------------------------------------------------------------------
if(${CMAKE_SYSTEM_NAME} MATCHES "Windows")
......@@ -42,11 +45,17 @@ endif(${CMAKE_SYSTEM_NAME} MATCHES "Windows")
option(BUILD_TESTING "To build the unit tests" ON)
option(BUILD_NR_SLICER_EXT "Where NiftyReg meets 3DSlicer" ${NiftyRegExtension_BUILD_SLICER_EXTENSION})
mark_as_advanced(BUILD_NR_SLICER_EXT)
if(NOT WIN32)
option(BUILD_ALL_DEP "All the dependencies are build" OFF)
else(NOT WIN32)
set(BUILD_ALL_DEP "TRUE")
endif(NOT WIN32)
option(USE_SSE "To enable SEE computation in some case" ON)
option(USE_DOUBLE "To force double precision instead of single precision" OFF)
option(USE_OPENMP "To use openMP for multi-CPU processing" ON)
option(BUILD_SHARED_LIBS "Build the libraries as shared" OFF)
option(USE_CUDA "To enable CUDA for a GPU implementation of the code" OFF)
option(BUILD_ALL_DEP "All the dependencies are build" OFF)
#-----------------------------------------------------------------------------
if(WIN32)
set(BUILD_ALL_DEP ON CACHE BOOL "All the dependencies are build" FORCE)
endif(WIN32)
#-----------------------------------------------------------------------------
# All dependencies are build to create the 3DSlicer package
if(BUILD_NR_SLICER_EXT)
set(BUILD_ALL_DEP ON)
......@@ -54,38 +63,6 @@ if(BUILD_NR_SLICER_EXT)
else(BUILD_NR_SLICER_EXT)
mark_as_advanced(CLEAR BUILD_ALL_DEP)
endif(BUILD_NR_SLICER_EXT)
option(USE_SSE "To enable SEE computation in some case" ON)
option(USE_DOUBLE "To force double precision instead of single precision" OFF)
option(USE_OPENMP "To use openMP for multi-CPU processing" ON)
if(USE_OPENMP)
find_package(OpenMP)
if (NOT OPENMP_FOUND)
set(USE_OPENMP OFF CACHE BOOL "To use openMP for multi-CPU processing" FORCE)
message(FATAL "OpenMP not found, forcing USE_OPENMP to OFF")
endif()
endif()
option(BUILD_SHARED_LIBS "Build the libraries as shared." OFF)
if(BUILD_SHARED_LIBS)
if(USE_CUDA)
set(BUILD_SHARED_LIBS OFF CACHE BOOL "Build the libraries as shared." FORCE)
set(NIFTYREG_LIBRARY_TYPE STATIC)
else()
set(NIFTYREG_LIBRARY_TYPE SHARED)
endif()
else()
set(NIFTYREG_LIBRARY_TYPE STATIC)
endif()
#-----------------------------------------------------------------------------
include_directories(${CMAKE_SOURCE_DIR}/reg-lib)
include_directories(${CMAKE_SOURCE_DIR}/reg-io)
include_directories(${CMAKE_SOURCE_DIR}/reg-io/nifti)
include_directories(${CMAKE_SOURCE_DIR}/third-party)
#-----------------------------------------------------------------------------
# Z library
# Try first to find the z library on the system and built is from the sources if it can not be find
......@@ -133,9 +110,11 @@ else(NOT BUILD_ALL_DEP)
endif(NOT BUILD_ALL_DEP)
include_directories(${CMAKE_SOURCE_DIR}/reg-io/png)
include_directories(${PNG_INCLUDE_DIR})
#-----------------------------------------------------------------------------
# NRRD file format support - The nrrd file format has been embedded into the project
include_directories(${CMAKE_SOURCE_DIR}/reg-lib)
include_directories(${CMAKE_SOURCE_DIR}/reg-io)
include_directories(${CMAKE_SOURCE_DIR}/reg-io/nifti)
include_directories(${CMAKE_SOURCE_DIR}/third-party)
include_directories(${CMAKE_BINARY_DIR})
include_directories(${CMAKE_SOURCE_DIR}/reg-io/nrrd)
include_directories(${CMAKE_SOURCE_DIR}/reg-io/nrrd/NrrdIO)
......@@ -153,18 +132,35 @@ if(USE_DOUBLE)
endif(USE_CUDA)
add_definitions(-D_USE_NR_DOUBLE)
endif(USE_DOUBLE)
#-----------------------------------------------------------------------------
if(USE_SSE)
set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -msse")
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -msse")
add_definitions(-D_USE_SSE)
endif(USE_SSE)
#-----------------------------------------------------------------------------
if(USE_OPENMP)
find_package(OpenMP)
if(NOT OPENMP_FOUND)
set(USE_OPENMP OFF CACHE BOOL "To use openMP for multi-CPU processing" FORCE)
message(WARNING "OpenMP does not appear to be supported by your compiler, forcing USE_OPENMP to OFF")
else(NOT OPENMP_FOUND)
set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} ${OpenMP_C_FLAGS}")
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} ${OpenMP_CXX_FLAGS}")
endif(NOT OPENMP_FOUND)
endif(USE_OPENMP)
#-----------------------------------------------------------------------------
if(BUILD_SHARED_LIBS)
if(USE_CUDA)
set(BUILD_SHARED_LIBS OFF CACHE BOOL "Build the libraries as shared." FORCE)
message(WARNING "CUDA is not compatible with shared libraries. Forcing BUILD_SHARED_LIBS to OFF")
set(NIFTYREG_LIBRARY_TYPE STATIC)
else(USE_CUDA)
set(NIFTYREG_LIBRARY_TYPE SHARED)
endif(USE_CUDA)
else(BUILD_SHARED_LIBS)
set(NIFTYREG_LIBRARY_TYPE STATIC)
endif(BUILD_SHARED_LIBS)
#-----------------------------------------------------------------------------
if(USE_CUDA)
# Check if the CUDA driver are available
......@@ -185,17 +181,14 @@ if(USE_CUDA)
)
# Check if the executable could not compile
if(NOT COMPILE_RESULT_VAR)
message(STATUS "The code to check the presence of a CUDA-enabled card failed.")
message(STATUS "The USE_CUDA flag has been turned OFF.")
set(USE_CUDA false)
message(WARNING "The code to check the presence of a CUDA-enabled card failed.")
message(WARNING "The USE_CUDA flag has been turned OFF.")
set(USE_CUDA OFF CACHE BOOL "To enable CUDA for a GPU implementation of the code" FORCE)
else(NOT COMPILE_RESULT_VAR)
# Check if the executable return failure
if(RUN_RESULT_VAR)
message(STATUS "No CUDA-enabled card has been detected.")
message(STATUS "Result message: ${RUN_RESULT_VAR}")
message(STATUS "Error message: ${RUN_OUTPUT_VAR}")
message(STATUS "The USE_CUDA flag has been turned OFF.")
set(USE_CUDA false)
message(WARNING "No CUDA-enabled card has been detected\nResult message: ${RUN_RESULT_VAR}\nError message: ${RUN_OUTPUT_VAR}\nThe USE_CUDA flag has been turned OFF.")
set(USE_CUDA OFF CACHE BOOL "To enable CUDA for a GPU implementation of the code" FORCE)
else(RUN_RESULT_VAR)
# Check if ptxas information should be displayed
set(CUDA_CAPABILITY "" CACHE STRING "CUDA capability to use to generate the ptxas information (1.0, 2.0, 3.0) [0]")
......@@ -216,7 +209,8 @@ if(USE_CUDA)
endif(RUN_RESULT_VAR)
endif(NOT COMPILE_RESULT_VAR)
else(CUDA_FOUND)
message(FATAL_ERROR "CUDA not found. Please turn the USE_CUDA flag off or set the variables manually.")
set(USE_CUDA OFF CACHE BOOL "To enable CUDA for a GPU implementation of the code" FORCE)
message(WARNING "CUDA is not compatible with shared libraries. Forcing BUILD_SHARED_LIBS to OFF")
endif(CUDA_FOUND)
endif(USE_CUDA)
#-----------------------------------------------------------------------------
......
......@@ -244,11 +244,21 @@ int main(int argc, char **argv)
}
else if(strcmp(argv[i], "-%v")==0 || strcmp(argv[i], "--vv")==0)
{
blockPercentage=atof(argv[++i]);
float value=atof(argv[++i]);
if(value<0.f || value>100.f){
reg_print_msg_error("The variance argument is expected to be between 0 and 100");
return EXIT_FAILURE;
}
blockPercentage=value;
}
else if(strcmp(argv[i], "-%i")==0 || strcmp(argv[i], "--ii")==0)
{
inlierLts=atof(argv[++i]);
float value=atof(argv[++i]);
if(value<0.f || value>100.f){
reg_print_msg_error("The inlier argument is expected to be between 0 and 100");
return EXIT_FAILURE;
}
inlierLts=value;
}
else if(strcmp(argv[i], "-interp")==0 || strcmp(argv[i], "--interp")==0)
{
......
......@@ -19,11 +19,12 @@ double reg_getKLDivergence1(nifti_image *referenceImage,
nifti_image *jacobianDetImg,
int *mask)
{
size_t voxelNumber = (size_t)referenceImage->nx*referenceImage->ny*referenceImage->nz;
#ifdef _WIN32
int voxel;
long voxel;
long voxelNumber = (long)referenceImage->nx*referenceImage->ny*referenceImage->nz;
#else
size_t voxel;
size_t voxel;
size_t voxelNumber = (size_t)referenceImage->nx*referenceImage->ny*referenceImage->nz;
#endif
DTYPE *refPtr=static_cast<DTYPE *>(referenceImage->data);
......@@ -138,11 +139,12 @@ void reg_getKLDivergenceVoxelBasedGradient1(nifti_image *referenceImage,
nifti_image *jacobianDetImg,
int *mask)
{
size_t voxelNumber = (size_t)referenceImage->nx*referenceImage->ny*referenceImage->nz;
#ifdef _WIN32
int voxel;
long voxel;
long voxelNumber = (long)referenceImage->nx*referenceImage->ny*referenceImage->nz;
#else
size_t voxel;
size_t voxelNumber = (size_t)referenceImage->nx*referenceImage->ny*referenceImage->nz;
#endif
DTYPE *refPtr=static_cast<DTYPE *>(referenceImage->data);
......
......@@ -84,8 +84,15 @@ double reg_getDTIMeasureValue(nifti_image *referenceImage,
float currentValue
)
{
#ifdef _WIN32
long voxel;
long voxelNumber = (long)referenceImage->nx*
referenceImage->ny*referenceImage->nz;
#else
size_t voxel;
size_t voxelNumber = (size_t)referenceImage->nx*
referenceImage->ny*referenceImage->nz;
#endif
/* As the tensor has 6 unique components that we need to worry about, read them out
for the floating and reference images. */
......@@ -115,13 +122,6 @@ double reg_getDTIMeasureValue(nifti_image *referenceImage,
// std::cerr << "BAM " << SUM_TEST << std::endl;
// reg_exit(1);
// Create some variables to be use in the openmp loop
#ifdef _WIN32
int voxel;
#else
size_t voxel;
#endif
double DTI_cost=0.0, n=0.0;
const double twoThirds = (2.0/3.0);
DTYPE rXX, rXY, rYY, rXZ, rYZ, rZZ;
......@@ -246,11 +246,12 @@ void reg_getVoxelBasedDTIMeasureGradient(nifti_image *referenceImage,
float currentValue)
{
// Create pointers to the reference and warped images
size_t voxelNumber = (size_t)referenceImage->nx*referenceImage->ny*referenceImage->nz;
#ifdef _WIN32
int voxel;
long voxel;
long voxelNumber = (long)referenceImage->nx*referenceImage->ny*referenceImage->nz;
#else
size_t voxel;
size_t voxel;
size_t voxelNumber = (size_t)referenceImage->nx*referenceImage->ny*referenceImage->nz;
#endif
/* As the tensor has 6 unique components that we need to worry about, read them out
......
......@@ -149,19 +149,19 @@ void reg_fem_getDeformationField(float *nodePositions,
float *femInterpolationWeight
)
{
size_t voxelNumber=
(size_t)deformationFieldImage->nx*
deformationFieldImage->ny*
deformationFieldImage->nz;
#ifdef _WIN32
long voxel;
long voxelNumber=(long)deformationFieldImage->nx*
deformationFieldImage->ny*deformationFieldImage->nz;
#else
size_t voxel;
size_t voxelNumber=(size_t)deformationFieldImage->nx*
deformationFieldImage->ny*deformationFieldImage->nz;
#endif
float *defPtrX = static_cast<float *>(deformationFieldImage->data);
float *defPtrY = &defPtrX[voxelNumber];
float *defPtrZ = &defPtrY[voxelNumber];
#ifdef _WIN32
int voxel;
#else
size_t voxel;
#endif
float coefficients[4];
float positionA[3], positionB[3], positionC[3], positionD[3];
#if defined (NDEBUG) && defined (_OPENMP)
......
......@@ -92,13 +92,21 @@ void reg_lncc::UpdateLocalStatImages(nifti_image *originalImage,
reg_tools_multiplyImageToImage(stdDevImage, stdDevImage, stdDevImage);
reg_tools_kernelConvolution(meanImage, this->kernelStandardDeviation, this->kernelType, mask, this->activeTimePoint);
reg_tools_kernelConvolution(stdDevImage, this->kernelStandardDeviation, this->kernelType, mask, this->activeTimePoint);
#ifdef _WIN32
long voxel;
long voxelNumber=(long)originalImage->nvox;
#else
size_t voxel;
size_t voxelNumber=originalImage->nvox;
#endif
#if defined (NDEBUG) && defined (_OPENMP)
#pragma omp parallel for default(none) \
shared(originalImage, sdevPtr, meanPtr) \
shared(voxelNumber, sdevPtr, meanPtr) \
private(voxel)
#endif
for(voxel=0; voxel<originalImage->nvox; ++voxel)
for(voxel=0; voxel<voxelNumber; ++voxel)
{
// G*(I^2) - (G*I)^2
sdevPtr[voxel] = sqrt(sdevPtr[voxel] - reg_pow2(meanPtr[voxel]));
......@@ -285,8 +293,15 @@ double reg_getLNCCValue(nifti_image *referenceImage,
DTYPE *refSdevPtr=static_cast<DTYPE *>(referenceSdevImage->data);
DTYPE *warSdevPtr=static_cast<DTYPE *>(warpedSdevImage->data);
DTYPE *correlaPtr=static_cast<DTYPE *>(correlationImage->data);
int voxel, voxelNumber = (int)referenceImage->nx*
referenceImage->ny*referenceImage->nz;
#ifdef _WIN32
long voxel;
long voxelNumber=(long)referenceImage->nx*
referenceImage->ny*referenceImage->nz;
#else
size_t voxel;
size_t voxelNumber=(size_t)referenceImage->nx*
referenceImage->ny*referenceImage->nz;
#endif
// Iteration over all time points
for(int t=0; t<referenceImage->nt; ++t)
......@@ -460,8 +475,15 @@ void reg_getVoxelBasedLNCCGradient(nifti_image *referenceImage,
DTYPE *warSdevPtr=static_cast<DTYPE *>(warpedSdevImage->data);
DTYPE *correlaPtr=static_cast<DTYPE *>(correlationImage->data);
size_t voxel, voxelNumber = (int)referenceImage->nx *
referenceImage->ny * referenceImage->nz;
#ifdef _WIN32
long voxel;
long voxelNumber=(long)referenceImage->nx*
referenceImage->ny*referenceImage->nz;
#else
size_t voxel;
size_t voxelNumber=(size_t)referenceImage->nx*
referenceImage->ny*referenceImage->nz;
#endif
// Create some pointers to the gradient images
DTYPE *lnccGradPtrX = static_cast<DTYPE *>(lnccGradientImage->data);
......@@ -576,12 +598,17 @@ void reg_getVoxelBasedLNCCGradient(nifti_image *referenceImage,
}
// Check for NaN
DTYPE val;
#ifdef _WIN32
voxelNumber=(long)lnccGradientImage->nvox;
#else
voxelNumber=lnccGradientImage->nvox;
#endif
#if defined (NDEBUG) && defined (_OPENMP)
#pragma omp parallel for default(none) \
shared(lnccGradientImage,lnccGradPtrX) \
shared(voxelNumber,lnccGradPtrX) \
private(voxel, val)
#endif
for(voxel=0; voxel<lnccGradientImage->nvox; ++voxel)
for(voxel=0; voxel<voxelNumber; ++voxel)
{
val=lnccGradPtrX[voxel];
if(val!=val || isinf(val)!=0)
......
......@@ -2459,7 +2459,13 @@ void reg_defField_compose2D(nifti_image *deformationField,
int *mask)
{
size_t DFVoxelNumber=(size_t)deformationField->nx*deformationField->ny;
#ifdef _WIN32
long i;
long warVoxelNumber=(size_t)dfToUpdate->nx*dfToUpdate->ny;
#else
size_t i;
size_t warVoxelNumber=(size_t)dfToUpdate->nx*dfToUpdate->ny;
#endif
DTYPE *defPtrX = static_cast<DTYPE *>(deformationField->data);
DTYPE *defPtrY = &defPtrX[DFVoxelNumber];
......@@ -2479,11 +2485,6 @@ void reg_defField_compose2D(nifti_image *deformationField,
df_voxel2Real=&(deformationField->qto_xyz);
}
#ifdef _WIN32
int i;
#else
size_t i;
#endif
size_t index;
int a, b, pre[2];
DTYPE realDefX, realDefY, voxelX, voxelY;
......@@ -2562,7 +2563,13 @@ void reg_defField_compose3D(nifti_image *deformationField,
{
const int DefFieldDim[3]= {deformationField->nx,deformationField->ny,deformationField->nz};
const size_t DFVoxelNumber=(size_t)DefFieldDim[0]*DefFieldDim[1]*DefFieldDim[2];
size_t warVoxelNumber=dfToUpdate->nx*dfToUpdate->ny*dfToUpdate->nz;
#ifdef _WIN32
long i;
long warVoxelNumber=(size_t)dfToUpdate->nx*dfToUpdate->ny;
#else
size_t i;
size_t warVoxelNumber=(size_t)dfToUpdate->nx*dfToUpdate->ny;
#endif
DTYPE *defPtrX = static_cast<DTYPE *>(deformationField->data);
DTYPE *defPtrY = &defPtrX[DFVoxelNumber];
......@@ -2589,11 +2596,6 @@ void reg_defField_compose3D(nifti_image *deformationField,
df_voxel2Real=&deformationField->qto_xyz;
}
#ifdef _WIN32
int i;
#else
size_t i;
#endif
size_t tempIndex, index;
int a, b, c, currentX, currentY, currentZ, pre[3];
DTYPE realDef[3], voxel[3], basis, tempBasis;
......@@ -4137,11 +4139,14 @@ void compute_lie_bracket(nifti_image *img1,
)
{
reg_exit(1); // to update
#ifdef _WIN32
long voxelNumber=(long)img1->nx*img1->ny*img1->nz;
#else
size_t voxNumber=(size_t)img1->nx*img1->ny*img1->nz;
#endif
// Lie bracket using Jacobian for testing
if(use_jac)
{
size_t voxNumber = img1->nx*img1->ny*img1->nz;
mat33 *jacImg1=(mat33 *)malloc(voxNumber*sizeof(mat33));
mat33 *jacImg2=(mat33 *)malloc(voxNumber*sizeof(mat33));
......@@ -4259,17 +4264,19 @@ void compute_lie_bracket(nifti_image *img1,
// Compute the lie bracket value using difference of composition
#ifdef _WIN32
int i;
long i;
voxelNumber=(long)res->nvox;
#else
size_t i;
voxNumber=res->nvox;
#endif
#if defined (NDEBUG) && defined (_OPENMP)
#pragma omp parallel for default(none) \
shared(res, resPtr, one_twoPtr, two_onePtr) \
shared(voxNumber, resPtr, one_twoPtr, two_onePtr) \
private(i)
#endif
for(i=0; i<res->nvox; ++i)
for(i=0; i<voxNumber; ++i)
resPtr[i]=two_onePtr[i]-one_twoPtr[i];
// Free the temporary nifti images
nifti_image_free(one_two);
......@@ -4286,9 +4293,11 @@ void compute_BCH_update1(nifti_image *img1, // current field
DTYPE *res=(DTYPE *)malloc(img1->nvox*sizeof(DTYPE));
#ifdef _WIN32
int i;
long i;
long voxelNumber=(long)img1->nvox;
#else
size_t i;
size_t voxelNumber=img1->nvox;
#endif
bool use_jac=false;
......@@ -4298,10 +4307,10 @@ void compute_BCH_update1(nifti_image *img1, // current field
DTYPE *img2Ptr=static_cast<DTYPE *>(img2->data);
#if defined (NDEBUG) && defined (_OPENMP)
#pragma omp parallel for default(none) \
shared(img1,img1Ptr,img2Ptr, res) \
shared(voxelNumber,img1Ptr,img2Ptr, res) \
private(i)
#endif
for(i=0; i<img1->nvox; ++i)
for(i=0; i<voxelNumber; ++i)
res[i] = img1Ptr[i] + img2Ptr[i];
if(type>0)
......@@ -4316,10 +4325,10 @@ void compute_BCH_update1(nifti_image *img1, // current field
DTYPE *lie_bracket_img2_img1Ptr=static_cast<DTYPE *>(lie_bracket_img2_img1->data);
#if defined (NDEBUG) && defined (_OPENMP)
#pragma omp parallel for default(none) \
shared(img1, res, lie_bracket_img2_img1Ptr) \
shared(voxelNumber, res, lie_bracket_img2_img1Ptr) \
private(i)
#endif
for(i=0; i<img1->nvox; ++i)
for(i=0; i<voxelNumber; ++i)
res[i] += 0.5 * lie_bracket_img2_img1Ptr[i];
if(type>1)
......@@ -4331,10 +4340,10 @@ void compute_BCH_update1(nifti_image *img1, // current field
DTYPE *lie_bracket_img2_lie1Ptr=static_cast<DTYPE *>(lie_bracket_img2_lie1->data);
#if defined (NDEBUG) && defined (_OPENMP)
#pragma omp parallel for default(none) \
shared(img1, res, lie_bracket_img2_lie1Ptr) \
shared(voxelNumber, res, lie_bracket_img2_lie1Ptr) \
private(i)
#endif
for(i=0; i<img1->nvox; ++i)
for(i=0; i<voxelNumber; ++i)
res[i] += lie_bracket_img2_lie1Ptr[i]/12.0;
if(type>2)
......@@ -4346,10 +4355,10 @@ void compute_BCH_update1(nifti_image *img1, // current field
DTYPE *lie_bracket_img1_lie1Ptr=static_cast<DTYPE *>(lie_bracket_img1_lie1->data);
#if defined (NDEBUG) && defined (_OPENMP)
#pragma omp parallel for default(none) \
shared(img1, res, lie_bracket_img1_lie1Ptr) \
shared(voxelNumber, res, lie_bracket_img1_lie1Ptr) \
private(i)
#endif
for(i=0; i<img1->nvox; ++i)
for(i=0; i<voxelNumber; ++i)
res[i] -= lie_bracket_img1_lie1Ptr[i]/12.0;
nifti_image_free(lie_bracket_img1_lie1);
......@@ -4362,10 +4371,10 @@ void compute_BCH_update1(nifti_image *img1, // current field
DTYPE *lie_bracket_img1_lie2Ptr=static_cast<DTYPE *>(lie_bracket_img1_lie2->data);
#if defined (NDEBUG) && defined (_OPENMP)
#pragma omp parallel for default(none) \
shared(img1, res, lie_bracket_img1_lie2Ptr) \
shared(voxelNumber, res, lie_bracket_img1_lie2Ptr) \
private(i)
#endif
for(i=0; i<img1->nvox; ++i)
for(i=0; i<voxelNumber; ++i)
res[i] -= lie_bracket_img1_lie2Ptr[i]/24.0;
nifti_image_free(lie_bracket_img1_lie2);
}// >3
......
......@@ -30,17 +30,23 @@ void svd(T ** in, size_t size_m, size_t size_n, T * w, T ** v)
reg_exit(1);
}
#ifdef _WIN32
long sm, sn, sn2;
long size__m=(long)size_m,size__n=(long)size_n;
#else
size_t sm, sn, sn2;
size_t size__m=size_m,size__n=size_n;
#endif
Eigen::MatrixXf m(size_m,size_n);
#if defined (NDEBUG) && defined (_OPENMP)
#pragma omp parallel for default(none) \
shared(in,m, size_m, size_n) \
shared(in,m, size__m, size__n) \
private(sm, sn)
#endif
for(sm=0; sm<size_m; sm++)
for(sm=0; sm<size__m; sm++)
{
for(sn=0; sn<size_n; sn++)
for(sn=0; sn<size__n; sn++)
{
m(sm,sn)=static_cast<T>(in[sm][sn]);
}
......@@ -50,17 +56,17 @@ void svd(T ** in, size_t size_m, size_t size_n, T * w, T ** v)
#if defined (NDEBUG) && defined (_OPENMP)
#pragma omp parallel for default(none) \
shared(in,svd,v,w, size_n,size_m) \
shared(in,svd,v,w, size__n,size__m) \
private(sn2, sn, sm)
#endif
for(sn=0; sn<size_n; sn++)
for(sn=0; sn<size__n; sn++)
{
w[sn]=svd.singularValues()(sn);
for(sn2=0; sn2<size_n; sn2++)
for(sn2=0; sn2<size__n; sn2++)
{
v[sn2][sn]=static_cast<T>(svd.matrixV()(sn2,sn));
}
for(sm=0; sm<size_m; sm++)
for(sm=0; sm<size__m; sm++)
{
in[sm][sn]=static_cast<T>(svd.matrixU()(sm,sn));
}
......
......@@ -656,7 +656,13 @@ void reg_getVoxelBasedNMIGradient3D(nifti_image *referenceImage,
)
{
//
#ifdef WIN32
long i;
long voxelNumber = (long)referenceImage->nx*referenceImage->ny*referenceImage->nz;
#else
size_t i;
size_t voxelNumber = (size_t)referenceImage->nx*referenceImage->ny*referenceImage->nz;
#endif
// Pointers to the image data
DTYPE *refImagePtr = static_cast<DTYPE *>(referenceImage->data);
DTYPE *warImagePtr = static_cast<DTYPE *>(warpedImage->data);
......@@ -683,7 +689,6 @@ void reg_getVoxelBasedNMIGradient3D(nifti_image *referenceImage,
size_t referenceOffset=referenceBinNumber[t]*floatingBinNumber[t];
size_t floatingOffset=referenceOffset+referenceBinNumber[t];
int r,w;
size_t i;
DTYPE refValue,warValue,gradX,gradY,gradZ;
double jointDeriv[3],refDeriv[3],warDeriv[3],commun,jointLog,refLog,warLog;
// Iterate over all voxel
......
......@@ -307,13 +307,16 @@ template <class T>
void reg_conjugateGradient<T>::UpdateGradientValues()
{
#ifdef _WIN32
int i;
#ifdef WIN32
long i;
long num = (long)this->dofNumber;