Mercurial > hg > minc-tools
changeset 2636:b0e5228b9c74
Merge remote-tracking branch 'develop/develop' into develop
author | Vladimir S. FONOV <vladimir.fonov@gmail.com> |
---|---|
date | Tue, 13 Mar 2012 11:12:31 -0400 |
parents | 356c5bda3c9d (current diff) 148cf7a20ca4 (diff) |
children | 0b9a5b816213 |
files | minc4itk/examples/itk_convert.cpp minc4itk/examples/itk_distance.cpp minc4itk/examples/itk_resample.cpp |
diffstat | 18 files changed, 1349 insertions(+), 812 deletions(-) [+] |
line wrap: on
line diff
--- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -3,7 +3,7 @@ # Andrew Janke - a.janke@gmail.com # Vladimir S. FONOV - vladimir.fonov@gmail.com -CMAKE_MINIMUM_REQUIRED(VERSION 2.8) +CMAKE_MINIMUM_REQUIRED(VERSION 2.6) INCLUDE(ExternalProject) PROJECT(minc2) @@ -12,6 +12,14 @@ ENABLE_TESTING() INCLUDE(CTest) +IF(CMAKE_VERSION VERSION_GREATER 2.8.2) +INCLUDE(BuildNETCDF) +INCLUDE(BuildHDF5) +SET(ENABLE_SUPERBUILD ON) +ELSE(CMAKE_VERSION VERSION_GREATER 2.8.2) +SET(ENABLE_SUPERBUILD OFF) +ENDIF(CMAKE_VERSION VERSION_GREATER 2.8.2) + #ADD_SUBDIRECTORY( doc ) # Packaging defines @@ -29,68 +37,85 @@ SET(PACKAGE_VERSION "${CPACK_PACKAGE_VERSION_MAJOR}.${CPACK_PACKAGE_VERSION_MINOR}.${CPACK_PACKAGE_VERSION_PATCH}") SET(PACKAGE_STRING "${PACKAGE_NAME} ${PACKAGE_VERSION}") -OPTION(BUILD_MINC2 "Support minc2 file format" ON) -OPTION(BUILD_TOOLS "Build minc tools (mincreshape,mincresample, etc)" ON) -OPTION(BUILD_CONVERTERS "Build minc conversion programs (mnc2nii, nii2mnc , dcm2mnc...)" ON) -OPTION(BUILD_EZMINC "Build C++ interface library EZminc" ON) +OPTION(BUILD_MINC2 "Support minc2 file format" ON) +OPTION(BUILD_TOOLS "Build minc tools (mincreshape,mincresample, etc)" ON) +OPTION(BUILD_CONVERTERS "Build minc conversion programs (mnc2nii, nii2mnc , dcm2mnc...)" ON) +OPTION(BUILD_EZMINC "Build C++ interface library EZminc" ON) +OPTION(BUILD_SHARED_LIBS "Build minc2 with shared libraries." OFF) + +IF(ENABLE_SUPERBUILD) OPTION(USE_SYSTEM_NETCDF "Use System NETCDF > 3.6.3" ON) +ELSE(ENABLE_SUPERBUILD) +SET(USE_SYSTEM_NETCDF ON) +ENDIF(ENABLE_SUPERBUILD) -IF(BUILD_EZMINC) -OPTION(BUILD_MINC4ITK "Build ITK interface" ON) -ENDIF(BUILD_EZMINC) +FIND_PACKAGE(ITK QUIET) + +IF(BUILD_EZMINC AND ITK_FOUND) + OPTION(BUILD_MINC4ITK "Build ITK interface" ON) +ENDIF(BUILD_EZMINC AND ITK_FOUND) IF(BUILD_MINC2) - OPTION(USE_SYSTEM_HDF5 "Use System HDF5 > 1.8 " ON) + IF(ITK_FOUND AND ITK_VERSION_MAJOR VERSION_EQUAL 4) # check if using ITK HDF5 library + + IF(ITKHDF5_LOADED) + message(WARNING "ITK compiled with HDF5 support!") + SET(USE_ITK_HDF5 ON) + ENDIF(ITKHDF5_LOADED) + ELSE(ITK_FOUND AND ITK_VERSION_MAJOR VERSION_EQUAL 4) # check if using ITK HDF5 library + + IF(ENABLE_SUPERBUILD) + OPTION(USE_SYSTEM_HDF5 "Use System HDF5 > 1.8" ON) + ELSE(ENABLE_SUPERBUILD) + SET(USE_SYSTEM_HDF5 ON) + ENDIF(ENABLE_SUPERBUILD) + + ENDIF(ITK_FOUND AND ITK_VERSION_MAJOR VERSION_EQUAL 4) # check if using ITK HDF5 library + ENDIF(BUILD_MINC2) +IF(BUILD_SHARED_LIBS) + SET(LIBRARY_TYPE SHARED) + SET(LIBRARY_INSTALL LIBRARY) +ELSE(BUILD_SHARED_LIBS) + SET(LIBRARY_TYPE STATIC) + SET(LIBRARY_INSTALL ARCHIVE) +ENDIF(BUILD_SHARED_LIBS) + # external packages IF(USE_SYSTEM_NETCDF) FIND_PACKAGE(NETCDF REQUIRED) ELSE(USE_SYSTEM_NETCDF) - ExternalProject_Add(NETCDF - SOURCE_DIR NETCDF - URL "ftp://ftp.unidata.ucar.edu/pub/netcdf/netcdf-4.0.1.tar.gz" - URL_MD5 "a251453c5477599f050fa4e593295186" - BUILD_IN_SOURCE 1 - INSTALL_DIR "${CMAKE_INSTALL_PREFIX}" - BUILD_COMMAND make - INSTALL_COMMAND make install - CONFIGURE_COMMAND ./configure --prefix=${CMAKE_INSTALL_PREFIX} --with-pic --disable-netcdf4 --disable-hdf4 --disable-dap --disable-shared --disable-cxx --disable-f77 --disable-f90 --disable-examples --enable-v2 --disable-docs - ) - - SET(NETCDF_LIBRARY ${CMAKE_INSTALL_PREFIX}/lib/libnetcdf.a ) - SET(NETCDF_INCLUDE_DIR ${CMAKE_INSTALL_PREFIX}/include ) + build_netcdf(${CMAKE_INSTALL_PREFIX}) ENDIF(USE_SYSTEM_NETCDF) IF(BUILD_MINC2) - IF(USE_SYSTEM_HDF5 ) -# SET(HDF5_FIND_COMPONENTS "C") -# SET(HDF5_USE_STATIC_LIBRARIES ON) - FIND_PACKAGE(HDF5 REQUIRED) - ELSE(USE_SYSTEM_HDF5) + IF(USE_ITK_HDF5) + SET(HDF5_INCLUDE_DIR ${ITKHDF5_INCLUDE_DIRS}) + SET(HDF5_LIBRARY ${ITKHDF5_LIBRARIES}) + SET(ZLIB_LIBRARIES ${ITKZLIB_LIBRARIES}) + + ELSE(USE_ITK_HDF5) + FIND_PACKAGE(ZLIB REQUIRED) - ExternalProject_Add(HDF5 - SOURCE_DIR HDF5 - URL "http://www.hdfgroup.org/ftp/HDF5/releases/hdf5-1.8.7/src/hdf5-1.8.7.tar.gz" - URL_MD5 "37711d4bcb72997e93d495f97c76c33a" - BUILD_IN_SOURCE 1 - INSTALL_DIR "${CMAKE_INSTALL_PREFIX}" - BUILD_COMMAND make - INSTALL_COMMAND make install - CONFIGURE_COMMAND ./configure --prefix=${CMAKE_INSTALL_PREFIX} --with-pic --disable-shared --disable-cxx --disable-f77 --disable-f90 --disable-examples --disable-hl --disable-docs - ) - - SET(HDF5_INCLUDE_DIR ${CMAKE_INSTALL_PREFIX}/include ) - SET(HDF5_LIBRARY ${CMAKE_INSTALL_PREFIX}/lib/libhdf5.a ) - - ENDIF(USE_SYSTEM_HDF5) + IF(USE_SYSTEM_HDF5 ) + FIND_PACKAGE(HDF5 REQUIRED) + ELSE(USE_SYSTEM_HDF5) + build_hdf5(${CMAKE_INSTALL_PREFIX}) + ENDIF(USE_SYSTEM_HDF5) + + ENDIF(USE_ITK_HDF5) ENDIF(BUILD_MINC2) +ADD_DEFINITIONS(-DHAVE_CONFIG_H) -ADD_DEFINITIONS(-DHAVE_CONFIG_H) -ADD_DEFINITIONS(-fPIC) +# add for building relocatable library +IF(UNIX) + SET(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -fPIC") + SET(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -fPIC") +ENDIF(UNIX) # aliases SET(VERSION "${PACKAGE_VERSION}") @@ -131,8 +156,7 @@ CHECK_INCLUDE_FILES(inttypes.h HAVE_INTTYPES_H) CHECK_INCLUDE_FILES(string.h HAVE_STRING_H) CHECK_INCLUDE_FILES(strings.h HAVE_STRINGS_H) -CHECK_INCLUDE_FILES(pwd.h HAVE_PWD_H) - +CHECK_INCLUDE_FILES(pwd.h HAVE_PWD_H) IF(ZLIB_FOUND) @@ -242,11 +266,11 @@ SET(MINC2_LIBRARY minc2) SET(MINC2_LIBRARIES ${MINC2_LIBRARY} ${HDF5_LIBRARY} ${NETCDF_LIBRARY} ${ZLIB_LIBRARIES} m ) - ADD_LIBRARY(minc2 STATIC ${minc_LIB_SRCS} ) + ADD_LIBRARY(minc2 ${LIBRARY_TYPE} ${minc_LIB_SRCS} ) TARGET_LINK_LIBRARIES(minc2 ${NETCDF_LIBRARY} ${HDF5_LIBRARY} ${ZLIB_LIBRARIES} m ) - INSTALL(TARGETS minc2 ARCHIVE DESTINATION lib) + INSTALL(TARGETS minc2 ${LIBRARY_INSTALL} DESTINATION lib) INSTALL(FILES libsrc/minc.h libsrc/ParseArgv.h libsrc/voxel_loop.h @@ -268,9 +292,9 @@ ELSE(BUILD_MINC2) SET(minc_LIB_SRCS ${minc1_LIB_SRCS} ) - ADD_LIBRARY(minc STATIC ${minc1_LIB_SRCS} ) + ADD_LIBRARY(minc ${LIBRARY_TYPE} ${minc1_LIB_SRCS} ) TARGET_LINK_LIBRARIES(minc ${NETCDF_LIBRARY} ) - INSTALL(TARGETS minc ARCHIVE DESTINATION lib) + INSTALL(TARGETS minc ${LIBRARY_INSTALL} DESTINATION lib) INSTALL(FILES libsrc/minc.h libsrc/ParseArgv.h libsrc/voxel_loop.h @@ -289,7 +313,6 @@ ADD_DEPENDENCIES(${MINC2_LIBRARY} NETCDF) ENDIF(NOT USE_SYSTEM_NETCDF) - # build the main minc2 library # volume_io2 @@ -332,15 +355,15 @@ IF(BUILD_MINC2) # build and install volume_io2 library - ADD_LIBRARY(volume_io2 STATIC ${volume_io_LIB_SRCS}) - INSTALL(TARGETS volume_io2 ARCHIVE DESTINATION lib) + ADD_LIBRARY(volume_io2 ${LIBRARY_TYPE} ${volume_io_LIB_SRCS}) + INSTALL(TARGETS volume_io2 ${LIBRARY_INSTALL} DESTINATION lib) SET_TARGET_PROPERTIES(volume_io2 PROPERTIES VERSION ${PACKAGE_VERSION} SOVERSION ${CPACK_PACKAGE_VERSION_MAJOR}) SET( VOLUME_IO_LIBRARY volume_io2) ADD_DEPENDENCIES(volume_io2 ${MINC2_LIBRARY}) ELSE(BUILD_MINC2) - ADD_LIBRARY(volume_io STATIC ${volume_io_LIB_SRCS}) - INSTALL(TARGETS volume_io ARCHIVE DESTINATION lib) + ADD_LIBRARY(volume_io ${LIBRARY_TYPE} ${volume_io_LIB_SRCS}) + INSTALL(TARGETS volume_io ${LIBRARY_INSTALL} DESTINATION lib) SET_TARGET_PROPERTIES(volume_io PROPERTIES VERSION ${PACKAGE_VERSION} SOVERSION ${CPACK_PACKAGE_VERSION_MAJOR}) SET( VOLUME_IO_LIBRARY volume_io) ADD_DEPENDENCIES(volume_io ${MINC2_LIBRARY}) @@ -383,9 +406,9 @@ ENDIF(BUILD_MINC4ITK AND BUILD_EZMINC) # config for the build directory -set(MINC2_USE_FILE_CONFIG ${CMAKE_CURRENT_BINARY_DIR}/UseMINC2.cmake) +set(MINC2_USE_FILE_CONFIG ${CMAKE_CURRENT_BINARY_DIR}/UseMINC2.cmake) -set(MINC2_INCLUDE_DIRS_CONFIG +set(MINC2_INCLUDE_DIRS_CONFIG ${CMAKE_CURRENT_SOURCE_DIR} ${CMAKE_SOURCE_DIR}/libsrc ${CMAKE_SOURCE_DIR}/volume_io/Include @@ -398,6 +421,14 @@ ) ENDIF(BUILD_MINC2) +IF(BUILD_EZMINC) + SET(EZMINC_LIBRARIES minc_io ${MINC2_LIBRARIES}) +ENDIF(BUILD_EZMINC) + +IF(BUILD_MINC4ITK) + SET(MINC4ITK_LIBRARIES minc4itk ${EZMINC_LIBRARIES}) +ENDIF(BUILD_MINC4ITK) + set(MINC2_LIBRARY_DIRS_CONFIG ${CMAKE_CURRENT_BINARY_DIR}) configure_file(MINC2Config.cmake.in @@ -409,7 +440,7 @@ # config for install dir set(MINC2_USE_FILE_CONFIG "${CMAKE_INSTALL_PREFIX}/lib/UseMINC2.cmake") -set(MINC2_INCLUDE_DIRS_CONFIG +set(MINC2_INCLUDE_DIRS_CONFIG ${CMAKE_INSTALL_PREFIX}/include ) @@ -421,7 +452,7 @@ configure_file(UseMINC2.cmake.in ${CMAKE_CURRENT_BINARY_DIR}/CMakeFiles/UseMINC2.cmake @ONLY) -install(FILES ${CMAKE_CURRENT_BINARY_DIR}/CMakeFiles/UseMINC2.cmake ${CMAKE_CURRENT_BINARY_DIR}/CMakeFiles/MINC2Config.cmake +INSTALL(FILES ${CMAKE_CURRENT_BINARY_DIR}/CMakeFiles/UseMINC2.cmake ${CMAKE_CURRENT_BINARY_DIR}/CMakeFiles/MINC2Config.cmake DESTINATION lib COMPONENT Development) @@ -430,6 +461,3 @@ IF(BUILD_TESTING) ADD_SUBDIRECTORY( testdir ) ENDIF(BUILD_TESTING) - - -
--- a/MINC2Config.cmake.in +++ b/MINC2Config.cmake.in @@ -1,15 +1,20 @@ set(HAVE_MINC2 @BUILD_MINC2@) set(HAVE_EZMINC @BUILD_EZMINC@) set(HAVE_MINC4ITK @BUILD_MINC4ITK@) - +set(USE_ITK_HDF5 @USE_ITK_HDF5@) set(NETCDF_INCLUDE_DIR "@NETCDF_INCLUDE_DIR@") set(HDF5_INCLUDE_DIR "@HDF5_INCLUDE_DIR@") set(NETCDF_LIBRARY "@NETCDF_LIBRARY@") set(HDF5_LIBRARY "@HDF5_LIBRARY@") set(ZLIB_LIBRARIES "@ZLIB_LIBRARIES@") +set(ITK_DIR "@ITK_DIR@") set(MINC2_INCLUDE_DIRS "@MINC2_INCLUDE_DIRS_CONFIG@") set(MINC2_LIBRARY_DIRS "@MINC2_LIBRARY_DIRS_CONFIG@") set(MINC2_USE_FILE "@MINC2_USE_FILE_CONFIG@") set(MINC2_LIBRARIES "@MINC2_LIBRARIES@") + +set(EZMINC_LIBRARIES "@EZMINC_LIBRARIES@") +set(MINC4ITK_LIBRARIES "@MINC4ITK_LIBRARIES@") +set(VOLUME_IO_LIBRARIES "@VOLUME_IO_LIBRARY@")
--- a/UseMINC2.cmake.in +++ b/UseMINC2.cmake.in @@ -1,21 +1,16 @@ -include_directories(${MINC2_INCLUDE_DIRS}) -link_directories(${MINC2_LIBRARY_DIRS}) - - -#FIND_PACKAGE(NETCDF REQUIRED) +INCLUDE_DIRECTORIES(${MINC2_INCLUDE_DIRS}) +LINK_DIRECTORIES(${MINC2_LIBRARY_DIRS}) INCLUDE_DIRECTORIES( ${NETCDF_INCLUDE_DIR} ) -if(HAVE_MINC2) - #FIND_PACKAGE(HDF5 REQUIRED) - #FIND_PACKAGE(ZLIB REQUIRED) - - INCLUDE_DIRECTORIES( ${HDF5_INCLUDE_DIR} ) - +IF(HAVE_MINC2) SET(MINC2 "1") ADD_DEFINITIONS( -DMINC2=1 ) - #SET(MINC_LIBRARIES volume_io2 minc2 netcdf hdf5 z) -else(HAVE_MINC2) - #SET(MINC_LIBRARIES volume_io minc netcdf ) -endif(HAVE_MINC2) + INCLUDE_DIRECTORIES(${MINC2_INCLUDE_DIRS}) + LINK_DIRECTORIES(${MINC2_LIBRARY_DIRS}) + INCLUDE_DIRECTORIES( ${HDF5_INCLUDE_DIR} ) +ENDIF(HAVE_MINC2) +IF(USE_ITK_HDF5) + LINK_DIRECTORIES(${ITK_LIBRARY_DIRS}) +ENDIF(USE_ITK_HDF5)
new file mode 100644 --- /dev/null +++ b/cmake-modules/BuildHDF5.cmake @@ -0,0 +1,18 @@ +macro(build_hdf5 install_prefix) + +ExternalProject_Add(HDF5 + SOURCE_DIR HDF5 + URL "http://www.hdfgroup.org/ftp/HDF5/releases/hdf5-1.8.7/src/hdf5-1.8.7.tar.gz" + URL_MD5 "37711d4bcb72997e93d495f97c76c33a" + BUILD_IN_SOURCE 1 + INSTALL_DIR "${install_prefix}" + BUILD_COMMAND make + INSTALL_COMMAND make install + CONFIGURE_COMMAND ./configure --prefix=${install_prefix} --with-pic --disable-shared --disable-cxx --disable-f77 --disable-f90 --disable-examples --disable-hl --disable-docs + ) + +SET(HDF5_INCLUDE_DIR ${install_prefix}/include ) +SET(HDF5_LIBRARY ${install_prefix}/lib/libhdf5.a ) + + +endmacro(build_hdf5) \ No newline at end of file
new file mode 100644 --- /dev/null +++ b/cmake-modules/BuildNETCDF.cmake @@ -0,0 +1,18 @@ + +macro(build_netcdf install_prefix) + +ExternalProject_Add(NETCDF + SOURCE_DIR NETCDF + URL "ftp://ftp.unidata.ucar.edu/pub/netcdf/netcdf-4.0.1.tar.gz" + URL_MD5 "a251453c5477599f050fa4e593295186" + BUILD_IN_SOURCE 1 + INSTALL_DIR "${install_prefix}" + BUILD_COMMAND make + INSTALL_COMMAND make install + CONFIGURE_COMMAND ./configure --prefix=${install_prefix} --with-pic --disable-netcdf4 --disable-hdf4 --disable-dap --disable-shared --disable-cxx --disable-f77 --disable-f90 --disable-examples --enable-v2 --disable-docs +) + +SET(NETCDF_LIBRARY ${install_prefix}/lib/libnetcdf.a ) +SET(NETCDF_INCLUDE_DIR ${install_prefix}/include ) + +endmacro(build_netcdf) \ No newline at end of file
--- a/config.h.cmake +++ b/config.h.cmake @@ -49,3 +49,6 @@ #cmakedefine01 HAVE_ZLIB #cmakedefine01 HAVE_STRINGS_H #cmakedefine01 HAVE_STRING_H + +#define H5Acreate_vers 2 +
--- a/ezminc/CMakeLists.txt +++ b/ezminc/CMakeLists.txt @@ -21,10 +21,10 @@ ) -ADD_LIBRARY( minc_io ${MINC_IO_HEADERS} ${MINC_IO_SRC}) +ADD_LIBRARY( minc_io ${LIBRARY_TYPE} ${MINC_IO_HEADERS} ${MINC_IO_SRC}) TARGET_LINK_LIBRARIES(minc_io ${MINC2_LIBRARIES}) -INSTALL(TARGETS minc_io ARCHIVE DESTINATION lib) +INSTALL(TARGETS minc_io ${LIBRARY_INSTALL} DESTINATION lib) INSTALL(FILES ${MINC_IO_HEADERS} DESTINATION include)
--- a/minc4itk/CMakeLists.txt +++ b/minc4itk/CMakeLists.txt @@ -1,4 +1,3 @@ -FIND_PACKAGE(ITK REQUIRED) INCLUDE(${ITK_USE_FILE}) IF(NOT ITK_LIBRARIES) @@ -9,7 +8,11 @@ MESSAGE( FATAL_ERROR "ITK include directories are not found!") ENDIF(NOT ITK_INCLUDE_DIRS ) -#OPTION(BUILD_ITK_PLUGIN "Build ITK plugin" OFF) +IF(BUILD_SHARED_LIBS AND NOT ITK_BUILD_SHARED) + message(WARNING "Building Shared library but ITK is linked statically!") +ENDIF(BUILD_SHARED_LIBS AND NOT ITK_BUILD_SHARED) + +OPTION(BUILD_ITK_PLUGIN "Build ITK plugin" OFF) OPTION(BUILD_MINC4ITK_EXAMPLES "Build minc4itk examples" ON) LINK_DIRECTORIES(${ITK_LIBRARY_DIRS}) @@ -44,15 +47,24 @@ z ) #VF ITK plugin is disabled for now -#IF(BUILD_ITK_PLUGIN) -#ADD_LIBRARY( minc4itk_plugin SHARED minc4itk_plugin.cxx) -#TARGET_LINK_LIBRARIES( minc4itk_plugin minc4itk) -#INSTALL(TARGETS minc4itk_plugin LIBRARY DESTINATION lib) -#ENDIF(BUILD_ITK_PLUGIN) +IF(BUILD_ITK_PLUGIN) +message(WARNING "ITK IO plugin might not work as expected!") +ADD_LIBRARY( minc4itk_plugin SHARED minc4itk_plugin.cxx) +TARGET_LINK_LIBRARIES( minc4itk_plugin minc4itk) +INSTALL(TARGETS minc4itk_plugin LIBRARY DESTINATION lib) +ENDIF(BUILD_ITK_PLUGIN) -INSTALL(TARGETS minc4itk ARCHIVE DESTINATION lib) +INSTALL(TARGETS minc4itk ${LIBRARY_INSTALL} DESTINATION lib) INSTALL(FILES ${MINC4ITK_HEADERS} DESTINATION include) IF(BUILD_MINC4ITK_EXAMPLES) add_subdirectory(examples) -ENDIF(BUILD_MINC4ITK_EXAMPLES) \ No newline at end of file +ENDIF(BUILD_MINC4ITK_EXAMPLES) + +IF(BUILD_TOOLS) +add_subdirectory(tools) +ENDIF(BUILD_TOOLS) + +IF(BUILD_CONVERTERS) +add_subdirectory(conversion) +ENDIF(BUILD_CONVERTERS)
new file mode 100644 --- /dev/null +++ b/minc4itk/conversion/itk_convert.cpp @@ -0,0 +1,765 @@ +#include <iostream> + +#include <itkImage.h> +#include <itkImageFileReader.h> +#include <itkImageFileWriter.h> +#include <itkImageIOFactory.h> +#include <itkImageIOBase.h> +#include <itkFlipImageFilter.h> +#include <itkMetaDataObject.h> +#include <itkMetaDataDictionary.h> +#include <itkCastImageFilter.h> + +#include "itkMincImageIOFactory.h" +#include "itkMincImageIO.h" + +#include <time_stamp.h> // for creating minc style history entry +#include <getopt.h> + +#include <vector> +#include <algorithm> +#include <iostream> +#include <iomanip> + +#include "minc_helpers.h" + +typedef std::vector<double> double_vector; + + +void show_usage (const char *name) +{ + std::cerr + << "Usage: "<<name<<" <input> <output> " << std::endl + << "--clobber clobber output files"<<std::endl + << "--verbose be verbose"<<std::endl + << "--show-meta show meta information (if present)"<<std::endl + << " Spatial transformation flags:"<<std::endl + << "--inv-x invert X axis"<<std::endl + << "--inv-y invert Y axis"<<std::endl + << "--inv-z invert Z axis"<<std::endl + << "--center set origin to the center of the image"<<std::endl + << " DWI related flags "<<std::endl + << "--dwi - assume that we are dealing with DWI scan"<<std::endl + << "--use-b-matrix - convert b-matrix (if present) into DWI gradient directions and b-value, implies DWI"<<std::endl + << "--minc-to-nrrd Convert minc style to nrrd style "<<std::endl + << "--nrrd-to-minc Convert nrrd style to minc style "<<std::endl + << " Data conversion flags: "<<std::endl + << "--char - cast data as signed char"<<std::endl + << "--byte - cast data as unsigned char"<<std::endl + << "--short - cast data as signed short"<<std::endl + << "--ushort - cast data as unsigned short"<<std::endl + << "--int - cast data as signed int"<<std::endl + << "--uint - cast data as unsigned int"<<std::endl + << "--float - cast data as float"<<std::endl + << "--double - cast data as double"<<std::endl + << " MINC file format related flags:" << std::endl + << "--mbyte - save pixel value using signed char"<<std::endl + << "--mshort - save pixel value using short"<<std::endl + << "--mint - save pixel value using int"<<std::endl + << "--mfloat - save pixel value using float"<<std::endl + << "--mdouble - save pixel value using double"<<std::endl + +; +} + +typedef itk::ImageIOBase IOBase; +typedef itk::SmartPointer<IOBase> IOBasePointer; + +class ImageConverterBase +{ +protected: + bool assume_dti; + bool use_b_matrix; + bool dwi_flip_z; + std::string history; + bool verbose; + bool inv_x; + bool inv_y; + bool inv_z; + bool center; + bool show_meta; + int minc_type; + bool minc_to_nrrd; + bool nrrd_to_minc; +public: + void setup(bool _verbose=false, + bool _assume_dti=false, + bool _use_b_matrix=false, + bool _dwi_flip_z=false, + bool _inv_x=false,bool _inv_y=false,bool _inv_z=false, + bool _center=false, + bool _show_meta=false, + const std::string _history="", + int _minc_type=-1, + bool _minc_to_nrrd=false, + bool _nrrd_to_minc=false + ) + { + assume_dti=_assume_dti; + use_b_matrix=_use_b_matrix; + dwi_flip_z=_dwi_flip_z; + history=_history; + verbose=_verbose; + inv_x=_inv_x; + inv_y=_inv_y; + inv_z=_inv_z; + center=_center; + show_meta=_show_meta; + minc_type=_minc_type; + minc_to_nrrd=_minc_to_nrrd; + nrrd_to_minc=_nrrd_to_minc; + } + + virtual void load_and_save_image(IOBase* base,const char *fname,itk::ImageIOBase::IOComponentType oct)=0; + + virtual ~ImageConverterBase() + {} +}; + +template<class TInputImage> class ImageConverter: public ImageConverterBase +{ +protected: + + typedef TInputImage InputImageType; + typedef typename InputImageType::PixelType InputImagePixelType; + +public: + + ~ImageConverter() + { + + } + + ImageConverter() + { + setup(); + } + + void convert_meta_minc_to_nrrd(itk::MetaDataDictionary& dict) + { + // let's try converting DTI-related meta-information + double_vector bvalues,direction_x,direction_y,direction_z; + + if( itk::ExposeMetaData<double_vector>( dict , "acquisition:bvalues",bvalues) && + itk::ExposeMetaData<double_vector>( dict , "acquisition:direction_x",direction_x) && + itk::ExposeMetaData<double_vector>( dict , "acquisition:direction_y",direction_y) && + itk::ExposeMetaData<double_vector>( dict , "acquisition:direction_z",direction_z)) + { + //We have got MINC-style DTI metadata + if(bvalues.size()!=direction_x.size() || + bvalues.size()!=direction_y.size() || + bvalues.size()!=direction_z.size() ) + { + std::cerr<<"WARNING: Different number of components of DTI directions"<<std::endl + <<" Skipping DTI metadata conversion!"<<std::endl; + } else { + double max_b_value=*std::max_element(bvalues.begin(),bvalues.end()); + if(verbose) + std::cout<<"Found maximum B-value:"<<max_b_value<<std::endl; + + itk::EncapsulateMetaData<std::string>( dict, + "modality", + "DWMRI"); + + //TODO: find out if we can get this info from minc + // using identity for now + std::vector< std::vector< double > > measurement_frame(3,std::vector< double >(3)); + measurement_frame[0][0]=1.0; + measurement_frame[1][1]=1.0; + measurement_frame[2][2]=1.0; + + + itk::EncapsulateMetaData< std::vector< std::vector< double > > >( dict, + "NRRD_measurement frame", + measurement_frame); + + + //TODO: deal with NRRD_thicknesses ? + + std::ostringstream ossKey; + ossKey<<std::setw(9) << std::setiosflags(std::ios::fixed) + << std::setprecision(6) << std::setiosflags(std::ios::right) << max_b_value; + + + itk::EncapsulateMetaData<std::string>( dict, + "NRRD_centerings[0]","cell"); + itk::EncapsulateMetaData<std::string>( dict, + "NRRD_centerings[1]","cell"); + itk::EncapsulateMetaData<std::string>( dict, + "NRRD_centerings[2]","cell"); + + itk::EncapsulateMetaData<std::string>( dict, + "NRRD_kinds[0]","space"); + itk::EncapsulateMetaData<std::string>( dict, + "NRRD_kinds[1]","space"); + itk::EncapsulateMetaData<std::string>( dict, + "NRRD_kinds[2]","space"); + + + itk::EncapsulateMetaData<std::string>( dict, + "DWMRI_b-value", + ossKey.str()); + + + size_t zero_b_value_cnt=0; + for(size_t i=0;i<bvalues.size();i++) + { + double direction[3]={direction_x[i],direction_y[i],direction_z[i]}; + double bval=bvalues[i]; + double b_scale=1.0; + + if(fabs(bval)<1e-3) + zero_b_value_cnt++; + else + { + b_scale=::sqrt(bval/max_b_value); + + direction[0]*=b_scale; + direction[1]*=b_scale; + direction[2]*=(dwi_flip_z&&!use_b_matrix)?-b_scale:b_scale; + } + + std::ostringstream ossKey; + ossKey << "DWMRI_gradient_" << std::setw(4) << std::setfill('0') << i; + + std::ostringstream ossMetaString; + ossMetaString << std::setw(9) << std::setiosflags(std::ios::fixed) + << std::setprecision(6) << std::setiosflags(std::ios::right) + << direction[0] + << " " + << std::setw(9) << std::setiosflags(std::ios::fixed) + << std::setprecision(6) << std::setiosflags(std::ios::right) + << direction[1] + << " " + << std::setw(9) << std::setiosflags(std::ios::fixed) + << std::setprecision(6) << std::setiosflags(std::ios::right) + << direction[2]; + + // std::cout<<ossKey.str()<<ossMetaString.str()<<std::endl; + itk::EncapsulateMetaData<std::string>( dict, + ossKey.str(), ossMetaString.str() ); + } + if(verbose) + std::cout<<"Found "<<zero_b_value_cnt<<" Zero b-value components"<<std::endl; + } + } + } + + void convert_meta_nrrd_to_minc(itk::MetaDataDictionary& dict) + { + //1. get b-value + std::string b_value_,gradient_value_; + double bval=0.0,vect3d[3]; + double_vector bvalues,direction_x,direction_y,direction_z; + + //remove B-matrix, it is of no use now + //if(dict.HasKey( "acquisition:b_matrix" )) + + + if(itk::ExposeMetaData<std::string>( dict,"DWMRI_b-value", b_value_)) + { + int i=0; + while(true) + { + bval=atof(b_value_.c_str()); + + std::ostringstream ossKey; + ossKey << "DWMRI_gradient_" << std::setw(4) << std::setfill('0') << i++; + if(itk::ExposeMetaData<std::string>( dict,ossKey.str(), gradient_value_)) + { + std::istringstream iss(gradient_value_); + iss >> vect3d[0] >> vect3d[1] >> vect3d[2]; + + //? normalize to unit vector and modulate bvalue + double b_scale=vect3d[0]*vect3d[0]+vect3d[1]*vect3d[1]+vect3d[2]*vect3d[2]; + double vmag=::sqrt(b_scale); + + if(vmag<0.1) //Zero + { + bvalues.push_back(0.0); + direction_x.push_back(0.0); + direction_y.push_back(0.0); + direction_z.push_back(0.0); + } else { + bvalues.push_back(bval*b_scale); + direction_x.push_back(vect3d[0]/vmag); + direction_y.push_back(vect3d[1]/vmag); + direction_z.push_back(((dwi_flip_z&&!use_b_matrix)?-vect3d[2]:vect3d[2])/vmag); + } + } else { + break; + } + } + + if(!bvalues.empty()) + { + itk::EncapsulateMetaData<double_vector>( dict , "acquisition:bvalues",bvalues); + itk::EncapsulateMetaData<double_vector>( dict , "acquisition:direction_x",direction_x); + itk::EncapsulateMetaData<double_vector>( dict , "acquisition:direction_y",direction_y); + itk::EncapsulateMetaData<double_vector>( dict , "acquisition:direction_z",direction_z); + itk::EncapsulateMetaData<double>( dict , "acquisition:b_value", bval); + } + } + } + + + + double decompose_b_matrix(const double_vector& b_matrix, + double_vector & bvalues, + double_vector & direction_x, + double_vector & direction_y, + double_vector & direction_z) + { + size_t elements=b_matrix.size(); + if(elements%6) //should be divisible by 6 + std::cerr<<"Number of elements in b_matrix is not divisible by 6 , probably an error!"<<std::endl; + elements/=6; + + bvalues.resize(elements); + direction_x.resize(elements); + direction_y.resize(elements); + direction_z.resize(elements); + + vnl_matrix_fixed<double, 3, 3> bMatrix; + double max_bval=0; + + for(size_t i=0;i<elements;i++) + { + //Code from DicomToNrrdConverter + // UNC comments: The principal eigenvector of the bmatrix is to be extracted as + // it's the gradient direction and trace of the matrix is the b-value + + // UNC comments: Fill out the 3x3 bmatrix with the 6 components read from the + // DICOM header. + bMatrix[0][0] = b_matrix[i*6+0]; + bMatrix[0][1] = b_matrix[i*6+1]; + bMatrix[0][2] = b_matrix[i*6+2]; + bMatrix[1][1] = b_matrix[i*6+3]; + bMatrix[1][2] = b_matrix[i*6+4]; + bMatrix[2][2] = b_matrix[i*6+5]; + bMatrix[1][0] = bMatrix[0][1]; + bMatrix[2][0] = bMatrix[0][2]; + bMatrix[2][1] = bMatrix[1][2]; + + double bvalue = bMatrix[0][0] + bMatrix[1][1] + bMatrix[2][2]; + + if(bvalue>max_bval) max_bval=bvalue; + + // UNC comments: The b-value si the trace of the bmatrix + // UNC comments: Even if the bmatrix is null, the svd decomposition set the 1st eigenvector + // to (1,0,0). So we force the gradient direction to 0 if the bvalue is null + if(bvalue == 0.0 ) //TODO: threshold? + { + direction_x[i] = 0.0; + direction_y[i] = 0.0; + direction_z[i] = 0.0; + bvalues[i] = 0.0; + } else { + // UNC comments: Computing the decomposition + vnl_svd<double> svd(bMatrix); + + // UNC comments: Extracting the principal eigenvector i.e. the gradient direction + + direction_x[i] = svd.U(0,0); + direction_y[i] = svd.U(1,0); + direction_z[i] = dwi_flip_z?svd.U(2,0):svd.U(2,0); + bvalues[i] = bvalue; + } + } + return max_bval; + } + + template<class TOutputImage> void _load_and_save_image(IOBase* base,const char *fname) + { + typename itk::ImageFileReader<TInputImage >::Pointer reader = itk::ImageFileReader<TInputImage >::New(); + typename itk::FlipImageFilter<TInputImage >::Pointer flip=itk::FlipImageFilter<TInputImage >::New(); + + + reader->SetImageIO(base); + reader->SetFileName(base->GetFileName()); + reader->Update(); + + typename TInputImage::Pointer img=reader->GetOutput(); + itk::MetaDataDictionary &thisDic=img->GetMetaDataDictionary(); + + if(show_meta) + { + + //let's write some meta information if there is any + for(itk::MetaDataDictionary::ConstIterator it=thisDic.Begin();it!=thisDic.End();++it) + { + itk::MetaDataObjectBase *bs=(*it).second; + itk::MetaDataObject<std::string> * str=dynamic_cast<itk::MetaDataObject<std::string> *>(bs); + if(str) + std::cout<<(*it).first.c_str()<<" = "<< str->GetMetaDataObjectValue().c_str()<<std::endl; + else + std::cout<<(*it).first.c_str()<<" type: "<< typeid(*bs).name()<<std::endl; + } + } + + if( thisDic.HasKey( "acquisition:b_matrix" ) && use_b_matrix ) + { + double_vector bmatrix; + + if( itk::ExposeMetaData<double_vector>( thisDic , "acquisition:b_matrix",bmatrix)) + { + double_vector bvalues,direction_x,direction_y,direction_z; + double bval=decompose_b_matrix(bmatrix,bvalues,direction_x,direction_y,direction_z); + + if(!bvalues.empty()) + { + itk::EncapsulateMetaData<double_vector>( thisDic , "acquisition:bvalues" ,bvalues); + itk::EncapsulateMetaData<double_vector>( thisDic , "acquisition:direction_x",direction_x); + itk::EncapsulateMetaData<double_vector>( thisDic , "acquisition:direction_y",direction_y); + itk::EncapsulateMetaData<double_vector>( thisDic , "acquisition:direction_z",direction_z); + itk::EncapsulateMetaData<double>( thisDic , "acquisition:b_value" ,bval); + } + } + } + + //making sure that all vectors contain the same number of parameters (just in case) + if( thisDic.HasKey( "acquisition:bvalues" ) && + thisDic.HasKey( "acquisition:direction_x") && + thisDic.HasKey( "acquisition:direction_y") && + thisDic.HasKey( "acquisition:direction_z") && + !nrrd_to_minc) + { + if(verbose) + std::cout<<"Converting MINC-style DWI to NRRD style"<<std::endl; + convert_meta_minc_to_nrrd(thisDic); + } else if ( thisDic.HasKey("DWMRI_b-value") && !minc_to_nrrd ) { + //We have got NRRD-style DTI metadata + if(verbose) + std::cout<<"Converting NRRD-style DWI to MINC style"<<std::endl; + convert_meta_nrrd_to_minc(thisDic); + } else if( nrrd_to_minc || minc_to_nrrd) { + std::cerr<<"ERROR: requested DWI headers are missing!"<<std::endl; + } + + if(verbose) + std::cout<<"Writing "<<fname<<"..."<<std::endl; + + if(inv_x||inv_y||inv_z) + { + typename itk::FlipImageFilter<TInputImage >::FlipAxesArrayType arr; + arr[0]=inv_x; + arr[1]=inv_y; + arr[2]=inv_z; + flip->SetFlipAxes(arr); + flip->SetInput(img); + flip->Update(); + img=flip->GetOutput(); + } + + if(center)//move origin to the center of the image + { + typename TInputImage::RegionType r=img->GetLargestPossibleRegion(); + //std::vector<double> corner[3]; + + typename TInputImage::IndexType idx; + typename TInputImage::PointType c; + + idx[0]=r.GetIndex()[0]+r.GetSize()[0]/2.0; + idx[1]=r.GetIndex()[1]+r.GetSize()[1]/2.0; + idx[2]=r.GetIndex()[2]+r.GetSize()[2]/2.0; + + img->TransformIndexToPhysicalPoint(idx,c); + + typename TInputImage::PointType org=img->GetOrigin(); + + org[0]-=c[0]; + org[1]-=c[1]; + org[2]-=c[2]; + + img->SetOrigin(org); + } + + if(!history.empty()) + minc::append_history(img,history); + + typename itk::CastImageFilter< TInputImage, TOutputImage >::Pointer cast=itk::CastImageFilter< TInputImage, TOutputImage >::New(); + + cast->SetInput(img); + + typename itk::ImageFileWriter< TOutputImage >::Pointer writer = itk::ImageFileWriter<TOutputImage >::New(); + writer->SetFileName(fname); + cast->Update(); + + cast->GetOutput()->SetMetaDataDictionary(thisDic); + + if(minc_type!=-1) + minc::set_minc_storage_type(cast->GetOutput(),(nc_type)minc_type,minc_type!=NC_BYTE); //store byte as unsigned only + + writer->SetInput( cast->GetOutput() ); + writer->Update(); + } + + + virtual void load_and_save_image(IOBase* io,const char *fname,itk::ImageIOBase::IOComponentType oct) + { + switch(oct) + { + case itk::ImageIOBase::UCHAR: + this->_load_and_save_image<itk::VectorImage<unsigned char, InputImageType::ImageDimension> >(io,fname); + break; + case itk::ImageIOBase::CHAR: + this->_load_and_save_image<itk::VectorImage<char, InputImageType::ImageDimension> >(io,fname); + break; + case itk::ImageIOBase::USHORT: + _load_and_save_image<itk::VectorImage<unsigned short, InputImageType::ImageDimension> >(io,fname); + break; + case itk::ImageIOBase::SHORT: + _load_and_save_image<itk::VectorImage<short, InputImageType::ImageDimension> >(io,fname); + break; + case itk::ImageIOBase::INT: + _load_and_save_image<itk::VectorImage<int, InputImageType::ImageDimension> >(io,fname); + break; + case itk::ImageIOBase::UINT: + _load_and_save_image<itk::VectorImage<unsigned int, InputImageType::ImageDimension> >(io,fname); + break; + case itk::ImageIOBase::FLOAT: + _load_and_save_image<itk::VectorImage<float, InputImageType::ImageDimension> >(io,fname); + break; + case itk::ImageIOBase::DOUBLE: + _load_and_save_image<itk::VectorImage<double, InputImageType::ImageDimension> >(io,fname); + break; + default: + itk::ExceptionObject("Unsupported component type"); + } + } +}; + +int main(int argc,char **argv) +{ + int verbose=0; + int clobber=0; + int show_meta=0; + int inv_x=0,inv_y=0,inv_z=0,center=0; + int assume_dti=0; + char *_history = time_stamp(argc, argv); + std::string history=_history; + int use_b_matrix=0; + int dwi_flip_z=0; + int minc_type=-1; + int nrrd_to_minc=0; + int minc_to_nrrd=0; + + int store_char=0,store_uchar=0,store_short=0,store_ushort=0,store_float=0,store_int=0,store_uint=0,store_double=0; + + free(_history); + + static struct option long_options[] = { + {"verbose", no_argument, &verbose, 1}, + {"quiet", no_argument, &verbose, 0}, + {"clobber", no_argument, &clobber, 1}, + {"inv-x", no_argument, &inv_x, 1}, + {"inv-y", no_argument, &inv_y, 1}, + {"inv-z", no_argument, &inv_z, 1}, + {"center", no_argument, ¢er, 1}, + {"show-meta", no_argument, &show_meta, 1}, + {"dti", no_argument, &assume_dti, 1}, + {"dwi", no_argument, &assume_dti, 1}, + {"use-b-matrix",no_argument, &use_b_matrix, 1}, + {"nrrd-to-minc",no_argument, &nrrd_to_minc, 1}, + {"minc-to-nrrd",no_argument, &minc_to_nrrd, 1}, + + {"float", no_argument, &store_float, 1}, + {"double", no_argument, &store_double, 1}, + {"byte", no_argument, &store_uchar, 1}, + {"char", no_argument, &store_char, 1}, + {"short", no_argument, &store_short, 1}, + {"ushort", no_argument, &store_ushort, 1}, + {"int", no_argument, &store_int, 1}, + {"uint", no_argument, &store_uint, 1}, + + {"mfloat", no_argument, &minc_type, NC_FLOAT}, + {"mdouble", no_argument, &minc_type, NC_DOUBLE}, + {"mbyte", no_argument, &minc_type, NC_BYTE}, + {"mshort", no_argument, &minc_type, NC_SHORT}, + {"mint", no_argument, &minc_type, NC_INT}, + + {0, 0, 0, 0} + }; + + int c; + for (;;) + { + /* getopt_long stores the option index here. */ + int option_index = 0; + + c = getopt_long (argc, argv, "", long_options, &option_index); + + /* Detect the end of the options. */ + if (c == -1) + break; + + switch (c) + { + case 0: + break; + case '?': + /* getopt_long already printed an error message. */ + default: + show_usage (argv[0]); + return 1; + } + } + + if((argc - optind)<2) + { + show_usage(argv[0]); + return 1; + } + std::string input=argv[optind]; + std::string output=argv[optind+1]; + + if (!clobber && !access (output.c_str(), F_OK)) + { + std::cerr << output.c_str () << " Exists!" << std::endl; + return 1; + } + if(nrrd_to_minc || minc_to_nrrd) + dwi_flip_z=1; + + try + { + itk::RegisterMincIO(); + + if(use_b_matrix) assume_dti=1; + + /* READING */ + if( verbose ) + std::cout<<"Reading "<<input.c_str()<<"..."<<std::endl; + + //try to figure out what we have got + IOBasePointer io = itk::ImageIOFactory::CreateImageIO(input.c_str(), itk::ImageIOFactory::ReadMode ); + + if(!io) + throw itk::ExceptionObject("Unsupported image file type"); + + io->SetFileName(input.c_str()); + io->ReadImageInformation(); + + size_t nd = io->GetNumberOfDimensions(); + size_t nc = io->GetNumberOfComponents(); + itk::ImageIOBase::IOComponentType ct = io->GetComponentType(); + itk::ImageIOBase::IOComponentType oct = ct; + + + if(store_char) + oct=itk::ImageIOBase::CHAR; + else if(store_uchar) + oct=itk::ImageIOBase::UCHAR; + else if(store_short) + oct=itk::ImageIOBase::SHORT; + else if(store_ushort) + oct=itk::ImageIOBase::USHORT; + else if(store_int) + oct=itk::ImageIOBase::INT; + else if(store_uint) + oct=itk::ImageIOBase::UINT; + else if(store_float) + oct=itk::ImageIOBase::FLOAT; + else if(store_double) + oct=itk::ImageIOBase::DOUBLE; + + std::string ct_s = io->GetComponentTypeAsString(ct); + std::string oct_s = io->GetComponentTypeAsString(oct); + + + if(verbose) + { + std::cout<<"dimensions:"<< nd << std::endl + <<"components:"<< nc << std::endl + <<"input type:"<< ct_s.c_str() <<std::endl + <<"output type:"<<oct_s.c_str()<<std::endl; + } + + + ImageConverterBase *converter=NULL; + + if(nd==3 || assume_dti) + { + if(verbose) std::cout<<"Writing 3D image..."<<std::endl; + switch(ct) + { + case itk::ImageIOBase::UCHAR : + converter=new ImageConverter<itk::VectorImage<unsigned char, 3> >(); + break; + case itk::ImageIOBase::CHAR : + converter=new ImageConverter<itk::VectorImage<char, 3> >(); + break; + case itk::ImageIOBase::USHORT : + converter=new ImageConverter<itk::VectorImage<unsigned short, 3> >(); + break; + case itk::ImageIOBase::SHORT : + converter=new ImageConverter<itk::VectorImage<short, 3> >(); + break; + case itk::ImageIOBase::INT : + converter=new ImageConverter<itk::VectorImage<int, 3> >(); + break; + case itk::ImageIOBase::UINT: + converter=new ImageConverter<itk::VectorImage<unsigned int, 3> >(); + break; + case itk::ImageIOBase::FLOAT : + converter=new ImageConverter<itk::VectorImage<float, 3> >(); + break; + case itk::ImageIOBase::DOUBLE: + converter=new ImageConverter<itk::VectorImage<double, 3> >(); + break; + default: + itk::ExceptionObject("Unsupported component type"); + } + } else if(nd==4 ) { //|| (nd==3 && nc>1) + if(verbose) std::cout<<"Writing 4D image..."<<std::endl; + switch(ct) + { + case itk::ImageIOBase::UCHAR: + converter=new ImageConverter<itk::VectorImage<unsigned char, 4> >(); + break; + case itk::ImageIOBase::CHAR: + converter=new ImageConverter<itk::VectorImage<char, 4> >(); + break; + case itk::ImageIOBase::USHORT: + converter=new ImageConverter<itk::VectorImage<unsigned short, 4> >(); + break; + case itk::ImageIOBase::SHORT: + converter=new ImageConverter<itk::VectorImage<short, 4> >(); + break; + case itk::ImageIOBase::INT: + converter=new ImageConverter<itk::VectorImage<int, 4> >(); + break; + case itk::ImageIOBase::UINT: + converter=new ImageConverter<itk::VectorImage<unsigned int, 4> >(); + break; + case itk::ImageIOBase::FLOAT: + converter=new ImageConverter<itk::VectorImage<float, 4> >(); + break; + case itk::ImageIOBase::DOUBLE: + converter=new ImageConverter<itk::VectorImage<double, 4> >(); + break; + default: + itk::ExceptionObject("Unsupported component type"); + } + } + else { + throw itk::ExceptionObject("Unsupported number of dimensions"); + } + + if(converter) + { + converter->setup(verbose,assume_dti,use_b_matrix,dwi_flip_z,inv_x,inv_y,inv_z,center,show_meta,history,minc_type,minc_to_nrrd,nrrd_to_minc); + converter->load_and_save_image(io,output.c_str(),oct); + + delete converter; + } + } + + catch( itk::ExceptionObject & err ) + { + std::cerr << "ExceptionObject caught !" << std::endl; + std::cerr << err << std::endl; + return 2; + } + return 0; +}
--- a/minc4itk/examples/CMakeLists.txt +++ b/minc4itk/examples/CMakeLists.txt @@ -1,15 +1,6 @@ LINK_LIBRARIES( minc4itk ) - -ADD_EXECUTABLE(itk_convert itk_convert.cpp) -ADD_EXECUTABLE(itk_distance itk_distance.cpp) ADD_EXECUTABLE(itk_dti itk_dti.cpp) -ADD_EXECUTABLE(itk_resample itk_resample.cpp) ADD_EXECUTABLE(volume_2_csv volume_2_csv.cpp) -INSTALL(TARGETS - itk_convert - itk_distance - itk_resample - DESTINATION bin)
deleted file mode 100644 --- a/minc4itk/examples/itk_convert.cpp +++ /dev/null @@ -1,309 +0,0 @@ -#include <iostream> - -#include <itkImage.h> -#include <itkImageFileReader.h> -#include <itkImageFileWriter.h> -#include <itkImageIOFactory.h> -#include <itkImageIOBase.h> -#include <itkFlipImageFilter.h> -#include <itkMetaDataObject.h> -#include <itkMetaDataDictionary.h> - -#include "itkMincImageIOFactory.h" -#include "itkMincImageIO.h" - -#include <getopt.h> - -void show_usage (const char *name) -{ - std::cerr - << "Usage: "<<name<<" <input> <output> " << std::endl - << "--clobber clobber output files"<<std::endl - << "--verbose be verbose"<<std::endl - << "--inv-x invert X axis"<<std::endl - << "--inv-y invert Y axis"<<std::endl - << "--inv-z invert Z axis"<<std::endl - << "--center set origin to the center of the image"<<std::endl - << "--show-meta show meta information (if present)"<<std::endl -; -} -typedef itk::ImageIOBase IOBase; -typedef itk::SmartPointer<IOBase> IOBasePointer; - -template<class ImageType> void load_and_save_image(IOBase* base, - const char *fname, - bool inv_x=false, - bool inv_y=false, - bool inv_z=false, - bool center=false, - bool verbose=false, - bool show_meta=false) -{ - typename itk::ImageFileReader<ImageType >::Pointer reader = itk::ImageFileReader<ImageType >::New(); - typename itk::FlipImageFilter<ImageType >::Pointer flip=itk::FlipImageFilter<ImageType >::New(); - - reader->SetImageIO(base); - reader->SetFileName(base->GetFileName()); - reader->Update(); - - typename ImageType::Pointer img=reader->GetOutput(); - - if(show_meta) - { - itk::MetaDataDictionary &thisDic=img->GetMetaDataDictionary(); - - //let's write some meta information if there is any - for(itk::MetaDataDictionary::ConstIterator it=thisDic.Begin();it!=thisDic.End();++it) - { - itk::MetaDataObjectBase *bs=(*it).second; - itk::MetaDataObject<std::string> * str=dynamic_cast<itk::MetaDataObject<std::string> *>(bs); - if(str) - std::cout<<(*it).first.c_str()<<" = "<< str->GetMetaDataObjectValue().c_str()<<std::endl; - else - std::cout<<(*it).first.c_str()<<" type: "<< typeid(*bs).name()<<std::endl; - } - } - - //TODO: convert metadata to and from minc format for DTI - /* WRITING */ - if(verbose) - std::cout<<"Writing "<<fname<<"..."<<std::endl; - - if(inv_x||inv_y||inv_z) - { - typename itk::FlipImageFilter<ImageType >::FlipAxesArrayType arr; - arr[0]=inv_x; - arr[1]=inv_y; - arr[2]=inv_z; - flip->SetFlipAxes(arr); - flip->SetInput(img); - flip->Update(); - img=flip->GetOutput(); - } - - if(center)//move origin to the center of the image - { - typename ImageType::RegionType r=img->GetLargestPossibleRegion(); - //std::vector<double> corner[3]; - - typename ImageType::IndexType idx; - typename ImageType::PointType c; - - idx[0]=r.GetIndex()[0]+r.GetSize()[0]/2.0; - idx[1]=r.GetIndex()[1]+r.GetSize()[1]/2.0; - idx[2]=r.GetIndex()[2]+r.GetSize()[2]/2.0; - - img->TransformIndexToPhysicalPoint(idx,c); - - typename ImageType::PointType org=img->GetOrigin(); - - org[0]-=c[0]; - org[1]-=c[1]; - org[2]-=c[2]; - - img->SetOrigin(org); - } - - typename itk::ImageFileWriter< ImageType >::Pointer writer = itk::ImageFileWriter<ImageType >::New(); - writer->SetFileName(fname); - writer->SetInput( img ); - writer->Update(); - -} - -int main(int argc,char **argv) -{ - int verbose=0; - int clobber=0; - int show_meta=0; - int inv_x=0,inv_y=0,inv_z=0,center=0; - //char *history = time_stamp(argc, argv); //maybe we should free it afterwards - - static struct option long_options[] = { - {"verbose", no_argument, &verbose, 1}, - {"quiet", no_argument, &verbose, 0}, - {"clobber", no_argument, &clobber, 1}, - {"inv-x", no_argument, &inv_x, 1}, - {"inv-y", no_argument, &inv_y, 1}, - {"inv-z", no_argument, &inv_z, 1}, - {"center", no_argument, ¢er, 1}, - {"show-meta", no_argument, &show_meta, 1}, - {0, 0, 0, 0} - }; - - int c; - for (;;) - { - /* getopt_long stores the option index here. */ - int option_index = 0; - - c = getopt_long (argc, argv, "", long_options, &option_index); - - /* Detect the end of the options. */ - if (c == -1) - break; - - switch (c) - { - case 0: - break; - case '?': - /* getopt_long already printed an error message. */ - default: - show_usage (argv[0]); - return 1; - } - } - - if((argc - optind)<2) - { - show_usage(argv[0]); - return 1; - } - std::string input=argv[optind]; - std::string output=argv[optind+1]; - - if (!clobber && !access (output.c_str(), F_OK)) - { - std::cerr << output.c_str () << " Exists!" << std::endl; - return 1; - } - - try - { - itk::RegisterMincIO(); - - /* READING */ - if(verbose) - std::cout<<"Reading "<<input.c_str()<<"..."<<std::endl; - - //try to figure out what we have got - IOBasePointer io = itk::ImageIOFactory::CreateImageIO(input.c_str(), itk::ImageIOFactory::ReadMode ); - - if(!io) - throw itk::ExceptionObject("Unsupported image file type"); - - io->SetFileName(input.c_str()); - io->ReadImageInformation(); - - size_t nd = io->GetNumberOfDimensions(); - size_t nc = io->GetNumberOfComponents(); - itk::ImageIOBase::IOComponentType ct = io->GetComponentType(); - std::string ct_s = io->GetComponentTypeAsString(ct); - - if(verbose) - { - std::cout<<"dimensions:"<<nd<<std::endl - <<"components:"<<nc<<std::endl - <<"type:"<<ct_s.c_str()<<std::endl; - } - - if(nd==3 && nc==1) - { - if(verbose) std::cout<<"Writing 3D image..."<<std::endl; - switch(ct) - { - case itk::ImageIOBase::UCHAR : - load_and_save_image<itk::Image<unsigned char, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::CHAR : - load_and_save_image<itk::Image<char, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::USHORT : - load_and_save_image<itk::Image<unsigned short, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::SHORT : - load_and_save_image<itk::Image<short, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::INT : - load_and_save_image<itk::Image<int, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::UINT: - load_and_save_image<itk::Image<unsigned int, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::FLOAT : - load_and_save_image<itk::Image<float, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::DOUBLE: - load_and_save_image<itk::Image<double, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - default: - itk::ExceptionObject("Unsupported component type"); - } - } else if((nd==4 && nc==1)) { - if(verbose) std::cout<<"Writing 4D image..."<<std::endl; - switch(ct) - { - case itk::ImageIOBase::UCHAR: - load_and_save_image<itk::Image<unsigned char, 4> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::CHAR: - load_and_save_image<itk::Image<char, 4> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::USHORT: - load_and_save_image<itk::Image<unsigned short, 4> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::SHORT: - load_and_save_image<itk::Image<short, 4> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::INT: - load_and_save_image<itk::Image<int, 4> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::UINT: - load_and_save_image<itk::Image<unsigned int, 4> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::FLOAT: - load_and_save_image<itk::Image<float, 4> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::DOUBLE: - load_and_save_image<itk::Image<double, 4> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - default: - itk::ExceptionObject("Unsupported component type"); - } - } else if((nd==3 && nc>1)) { - if(verbose) std::cout<<"Writing multicomponent 3D image..."<<std::endl; - - switch(ct) - { - case itk::ImageIOBase::UCHAR: - load_and_save_image<itk::VectorImage<unsigned char, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::CHAR: - load_and_save_image<itk::VectorImage<char, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::USHORT: - load_and_save_image<itk::VectorImage<unsigned short, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::SHORT: - load_and_save_image<itk::VectorImage<short, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::INT: - load_and_save_image<itk::VectorImage<int, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::UINT: - load_and_save_image<itk::VectorImage<unsigned int, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::FLOAT: - load_and_save_image<itk::VectorImage<float, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - case itk::ImageIOBase::DOUBLE: - load_and_save_image<itk::VectorImage<double, 3> >(io,output.c_str(),inv_x,inv_y,inv_z,center,verbose,show_meta); - break; - default: - itk::ExceptionObject("Unsupported component type"); - } - } else { - throw itk::ExceptionObject("Unsupported number of dimensions"); - } - - } - - catch( itk::ExceptionObject & err ) - { - std::cerr << "ExceptionObject caught !" << std::endl; - std::cerr << err << std::endl; - return 2; - } - return 0; -}
deleted file mode 100644 --- a/minc4itk/examples/itk_distance.cpp +++ /dev/null @@ -1,139 +0,0 @@ -/* ----------------------------- MNI Header ----------------------------------- -@NAME : -@DESCRIPTION: -@COPYRIGHT : - Copyright 2006 Vladimir Fonov, McConnell Brain Imaging Centre, - Montreal Neurological Institute, McGill University. - Permission to use, copy, modify, and distribute this - software and its documentation for any purpose and without - fee is hereby granted, provided that the above copyright - notice appear in all copies. The author and McGill University - make no representations about the suitability of this - software for any purpose. It is provided "as is" without - express or implied warranty. ----------------------------------------------------------------------------- */ -#include <unistd.h> -#include <getopt.h> -#include <iostream> -#include <itkDanielssonDistanceMapImageFilter.h> -#include <itkSignedDanielssonDistanceMapImageFilter.h> -#include <time_stamp.h> // for creating minc style history entry -#include "itkMincImageIOFactory.h" -#include "itkMincImageIO.h" -#include "minc_helpers.h" - -//#include <minc_wrappers.h> - -using namespace minc; -using namespace std; -void show_usage(const char *name) -{ - std::cerr - << "Usage: "<<name<<" <input> <output> " << endl - << "--verbose be verbose " << endl - << "--clobber clobber output files"<<endl - << "--signed produce signed distance map"<<endl; -} - -int main (int argc, char **argv) -{ - - int verbose=1; - double sigma=0.5; - double keep=1.0; - int order=5; - int approx=0; - int ss=0; - int clobber=0; - char *history = time_stamp(argc, argv); - - //int voxel_neibourhood=5; - static struct option long_options[] = { - {"clobber", no_argument, &clobber, 1}, - {"verbose", no_argument, &verbose, 1}, - {"quiet", no_argument, &verbose, 0}, - {"signed",no_argument, &ss, 1}, - {0, 0, 0, 0} - }; - - for (;;) { - /* getopt_long stores the option index here. */ - int option_index = 0; - - int c = getopt_long (argc, argv, "vq", long_options, &option_index); - - /* Detect the end of the options. */ - if (c == -1) break; - - switch (c) - { - // case 'n': - // voxel_neibourhood=atoi(optarg);break; - case 0: - break; - case '?': - /* getopt_long already printed an error message. */ - default: - show_usage (argv[0]); - return 1; - } - } - if ((argc - optind) < 2) { - show_usage (argv[0]); - return 1; - } - std::string input_f=argv[optind], out_f=argv[optind+1]; - - // check if the file already present - if (!clobber && !access (out_f.c_str (), F_OK)) - { - cerr << out_f.c_str () << " Exists!" << endl; - return 1; - } - - try - { - itk::ObjectFactoryBase::RegisterFactory(itk::MincImageIOFactory::New()); - itk::ImageFileReader<minc::mask3d >::Pointer reader = itk::ImageFileReader<minc::mask3d >::New(); - - //initializing the reader - reader->SetFileName(input_f.c_str()); - reader->Update(); - - minc::mask3d::Pointer input=reader->GetOutput(); - minc::image3d::Pointer output; - - typedef itk::DanielssonDistanceMapImageFilter< minc::mask3d, minc::image3d > - DistanceMapFilter; - typedef itk::SignedDanielssonDistanceMapImageFilter< minc::mask3d, minc::image3d > - SignedDistanceMapFilter; - - if(ss) - { - SignedDistanceMapFilter::Pointer dist(SignedDistanceMapFilter::New()); - dist->SetInput(input); - dist->Update(); - output=dist->GetOutput(); - } else { - DistanceMapFilter::Pointer dist(DistanceMapFilter::New()); - dist->SetInput(input); - dist->Update(); - output=dist->GetOutput(); - } - - minc::copy_metadata(output,input); - minc::append_history(output,history); - free(history); - - itk::ImageFileWriter< minc::image3d >::Pointer writer = itk::ImageFileWriter<minc::image3d >::New(); - writer->SetFileName(out_f.c_str()); - writer->SetInput( output ); - writer->Update(); - - } catch (const minc::generic_error & err) { - cerr << "Got an error at:" << err.file () << ":" << err.line () << endl; - cerr << err.msg()<<endl; - return 1; - } - return 0; -}
--- a/minc4itk/examples/itk_dti.cpp +++ b/minc4itk/examples/itk_dti.cpp @@ -75,6 +75,8 @@ try { + itk::RegisterMincIO(); + std::cout<<"Reading "<<argv[1]<<"..."<<std::endl; typedef itk::MincImageIO ImageIOType; @@ -88,7 +90,7 @@ double_vector bvalues,direction_x,direction_y,direction_z; //making sure that all vcrtors contain the same number of parameters (just in case) - if(!itk::ExposeMetaData<double_vector>( img->GetMetaDataDictionary() , "acquisition:bvalues",bvalues) || + if( !itk::ExposeMetaData<double_vector>( img->GetMetaDataDictionary() , "acquisition:bvalues",bvalues) || !itk::ExposeMetaData<double_vector>( img->GetMetaDataDictionary() , "acquisition:direction_x",direction_x) || !itk::ExposeMetaData<double_vector>( img->GetMetaDataDictionary() , "acquisition:direction_y",direction_y) || !itk::ExposeMetaData<double_vector>( img->GetMetaDataDictionary() , "acquisition:direction_z",direction_z)) @@ -96,6 +98,7 @@ std::cerr<<"Image doesn't have information on DTI gradients, can't process!"<<std::endl; return 2; } + if(bvalues.size()!=direction_x.size() || bvalues.size()!=direction_y.size() || bvalues.size()!=direction_z.size() )
deleted file mode 100644 --- a/minc4itk/examples/itk_resample.cpp +++ /dev/null @@ -1,259 +0,0 @@ -/* ----------------------------- MNI Header ----------------------------------- -@NAME : itk_resample -@DESCRIPTION: an example of using spline itnerpolation with MINC xfm transform -@COPYRIGHT : - Copyright 2011 Vladimir Fonov, McConnell Brain Imaging Centre, - Montreal Neurological Institute, McGill University. - Permission to use, copy, modify, and distribute this - software and its documentation for any purpose and without - fee is hereby granted, provided that the above copyright - notice appear in all copies. The author and McGill University - make no representations about the suitability of this - software for any purpose. It is provided "as is" without - express or implied warranty. ----------------------------------------------------------------------------- */ -#include <iostream> -#include <fstream> -#include <getopt.h> -#include <vector> -#include <valarray> -#include <math.h> -#include <limits> -#include <unistd.h> - -#include <itkResampleImageFilter.h> -#include <itkAffineTransform.h> -#include <itkNearestNeighborInterpolateImageFunction.h> -#include <itkBSplineInterpolateImageFunction.h> -#include <minc_general_transform.h> - -#include <unistd.h> -#include <getopt.h> -#include <time_stamp.h> // for creating minc style history entry - - -#include <itkImageFileReader.h> -#include <itkImageFileWriter.h> -#include <itkImageIOFactory.h> - -#include "itkMincImageIOFactory.h" -#include "itkMincImageIO.h" -#include "minc_helpers.h" - -//typedef itk::MincImageIO ImageIOType; -typedef itk::BSplineInterpolateImageFunction< minc::image3d, double, double > InterpolatorType; -typedef itk::ResampleImageFilter<minc::image3d, minc::image3d> FilterType; -typedef minc::XfmTransform<double,3,3> TransformType; - -using namespace std; - -void show_usage (const char * prog) -{ - std::cerr - << "Usage: "<<prog<<" <input> <output.mnc> " << std::endl - << "--clobber overwrite files" << std::endl - << "--like <example> (default behaviour analogous to use_input_sampling)"<<std::endl - << "--transform <xfm_transform> "<<std::endl - << "--order <n> spline order, default 2 "<<std::endl - << "--uniformize <step> - will make a volume with uniform step size and no direction cosines" << std::endl - << "--invert_transform - apply inverted transform"<<std::endl; -} - -void generate_uniform_sampling(FilterType* flt, const minc::image3d* img,double step) -{ - //obtain physical coordinats of all image corners - minc::image3d::RegionType r=img->GetLargestPossibleRegion(); - std::vector<double> corner[3]; - for(int i=0;i<8;i++) - { - minc::image3d::IndexType idx; - minc::image3d::PointType c; - idx[0]=r.GetIndex()[0]+r.GetSize()[0]*(i%2); - idx[1]=r.GetIndex()[1]+r.GetSize()[1]*((i/2)%2); - idx[2]=r.GetIndex()[2]+r.GetSize()[2]*((i/4)%2); - img->TransformIndexToPhysicalPoint(idx,c); - for(int j=0;j<3;j++) - corner[j].push_back(c[j]); - } - minc::image3d::IndexType start; - FilterType::SizeType size; - FilterType::OriginPointType org; - minc::image3d::SpacingType spc; - spc.Fill(step); - for(int j=0;j<3;j++) - { - std::sort(corner[j].begin(),corner[j].end()); - size[j]=ceil((corner[j][7]-corner[j][0])/step); - org[j]=corner[j][0]; - } - minc::image3d::DirectionType identity; - identity.SetIdentity(); - flt->SetOutputDirection(identity); - start.Fill(0); - flt->SetOutputStartIndex(start); - flt->SetSize(size); - flt->SetOutputOrigin(org); - flt->SetOutputSpacing(spc); -} - -int main (int argc, char **argv) -{ - int verbose=0, clobber=0,skip_grid=0; - int order=2; - std::string like_f,xfm_f,output_f,input_f; - double uniformize=0.0; - int invert=0; - char *history = time_stamp(argc, argv); - - static struct option long_options[] = { - {"verbose", no_argument, &verbose, 1}, - {"quiet", no_argument, &verbose, 0}, - {"clobber", no_argument, &clobber, 1}, - {"like", required_argument, 0, 'l'}, - {"transform", required_argument, 0, 't'}, - {"order", required_argument, 0, 'o'}, - {"uniformize", required_argument, 0, 'u'}, - {"invert_transform", no_argument, &invert, 1}, - {0, 0, 0, 0} - }; - - for (;;) { - /* getopt_long stores the option index here. */ - int option_index = 0; - - int c = getopt_long (argc, argv, "vqcl:t:o:u:", long_options, &option_index); - - /* Detect the end of the options. */ - if (c == -1) break; - - switch (c) - { - case 0: - break; - case 'v': - cout << "Version: 0.1" << endl; - return 0; - case 'l': - like_f=optarg; break; - case 't': - xfm_f=optarg; break; - case 'o': - order=atoi(optarg);break; - case 'u': - uniformize=atof(optarg);break; - case '?': - /* getopt_long already printed an error message. */ - default: - show_usage (argv[0]); - return 1; - } - } - - if ((argc - optind) < 2) { - show_usage(argv[0]); - return 1; - } - input_f=argv[optind]; - output_f=argv[optind+1]; - - if (!clobber && !access (output_f.c_str (), F_OK)) - { - std::cerr << output_f.c_str () << " Exists!" << std::endl; - return 1; - } - - try - { - itk::ObjectFactoryBase::RegisterFactory(itk::MincImageIOFactory::New()); - itk::ImageFileReader<minc::image3d >::Pointer reader = itk::ImageFileReader<minc::image3d >::New(); - - //initializing the reader - reader->SetFileName(input_f.c_str()); - reader->Update(); - - minc::image3d::Pointer in=reader->GetOutput(); - - FilterType::Pointer filter = FilterType::New(); - - //creating coordinate transformation objects - TransformType::Pointer transform = TransformType::New(); - if(!xfm_f.empty()) - { - //reading a minc style xfm file - transform->OpenXfm(xfm_f.c_str()); - if(!invert) transform->Invert(); //should be inverted by default to walk through target space - filter->SetTransform( transform ); - } - - //creating the interpolator - InterpolatorType::Pointer interpolator = InterpolatorType::New(); - interpolator->SetSplineOrder(order); - filter->SetInterpolator( interpolator ); - filter->SetDefaultPixelValue( 0 ); - - //this is for processing using batch system - filter->SetNumberOfThreads(1); - - minc::image3d::Pointer like=0; - if(!like_f.empty()) - { - itk::ImageFileReader<minc::image3d >::Pointer reader = itk::ImageFileReader<minc::image3d >::New(); - reader->SetFileName(like_f.c_str()); - reader->Update(); - if(uniformize!=0.0) - { - generate_uniform_sampling(filter,reader->GetOutput(),uniformize); - } else { - filter->SetOutputParametersFromImage(reader->GetOutput()); - filter->SetOutputDirection(reader->GetOutput()->GetDirection()); - } - like=reader->GetOutput(); - like->DisconnectPipeline(); - } - else - { - if(uniformize!=0.0) - { - generate_uniform_sampling(filter,in,uniformize); - } else { - //we are using original sampling - filter->SetOutputParametersFromImage(in); - filter->SetOutputDirection(in->GetDirection()); - } - } - - filter->SetInput(in); - filter->Update(); - //copy the metadate information, for some reason it is not preserved - //filter->GetOutput()->SetMetaDataDictionary(reader->GetOutput()->GetMetaDataDictionary()); - minc::image3d::Pointer out=filter->GetOutput(); - minc::copy_metadata(out,in); - minc::append_history(out,history); - free(history); - - //correct dimension order - if(like.IsNotNull()) - minc::copy_dimorder(out,like); - - //generic file writer - itk::ImageFileWriter< minc::image3d >::Pointer writer = itk::ImageFileWriter<minc::image3d >::New(); - writer->SetFileName(output_f.c_str()); - - writer->SetInput( out ); - //writer->UseInputMetaDataDictionaryOn(); - - writer->Update(); - - return 0; - } catch (const minc::generic_error & err) { - cerr << "Got an error at:" << err.file () << ":" << err.line () << endl; - return 1; - } - catch( itk::ExceptionObject & err ) - { - std::cerr << "ExceptionObject caught !" << std::endl; - std::cerr << err << std::endl; - return 2; - } - return 0; -};
--- a/minc4itk/examples/volume_2_csv.cpp +++ b/minc4itk/examples/volume_2_csv.cpp @@ -93,6 +93,8 @@ try { + itk::RegisterMincIO(); + minc::mask3d::Pointer mask(minc::mask3d::New()); //creating a minc reader
new file mode 100644 --- /dev/null +++ b/minc4itk/tools/itk_distance.cpp @@ -0,0 +1,140 @@ +/* ----------------------------- MNI Header ----------------------------------- +@NAME : +@DESCRIPTION: +@COPYRIGHT : + Copyright 2006 Vladimir Fonov, McConnell Brain Imaging Centre, + Montreal Neurological Institute, McGill University. + Permission to use, copy, modify, and distribute this + software and its documentation for any purpose and without + fee is hereby granted, provided that the above copyright + notice appear in all copies. The author and McGill University + make no representations about the suitability of this + software for any purpose. It is provided "as is" without + express or implied warranty. +---------------------------------------------------------------------------- */ +#include <unistd.h> +#include <getopt.h> +#include <iostream> +#include <itkDanielssonDistanceMapImageFilter.h> +#include <itkSignedDanielssonDistanceMapImageFilter.h> +#include <time_stamp.h> // for creating minc style history entry +#include "itkMincImageIOFactory.h" +#include "itkMincImageIO.h" +#include "minc_helpers.h" + +//#include <minc_wrappers.h> + +using namespace minc; +using namespace std; +void show_usage(const char *name) +{ + std::cerr + << "Usage: "<<name<<" <input> <output> " << endl + << "--verbose be verbose " << endl + << "--clobber clobber output files"<<endl + << "--signed produce signed distance map"<<endl; +} + +int main (int argc, char **argv) +{ + + int verbose=1; + double sigma=0.5; + double keep=1.0; + int order=5; + int approx=0; + int ss=0; + int clobber=0; + char *history = time_stamp(argc, argv); + + //int voxel_neibourhood=5; + static struct option long_options[] = { + {"clobber", no_argument, &clobber, 1}, + {"verbose", no_argument, &verbose, 1}, + {"quiet", no_argument, &verbose, 0}, + {"signed",no_argument, &ss, 1}, + {0, 0, 0, 0} + }; + + for (;;) { + /* getopt_long stores the option index here. */ + int option_index = 0; + + int c = getopt_long (argc, argv, "vq", long_options, &option_index); + + /* Detect the end of the options. */ + if (c == -1) break; + + switch (c) + { + // case 'n': + // voxel_neibourhood=atoi(optarg);break; + case 0: + break; + case '?': + /* getopt_long already printed an error message. */ + default: + show_usage (argv[0]); + return 1; + } + } + if ((argc - optind) < 2) { + show_usage (argv[0]); + return 1; + } + std::string input_f=argv[optind], out_f=argv[optind+1]; + + // check if the file already present + if (!clobber && !access (out_f.c_str (), F_OK)) + { + cerr << out_f.c_str () << " Exists!" << endl; + return 1; + } + + try + { + itk::RegisterMincIO(); + + itk::ImageFileReader<minc::mask3d >::Pointer reader = itk::ImageFileReader<minc::mask3d >::New(); + + //initializing the reader + reader->SetFileName(input_f.c_str()); + reader->Update(); + + minc::mask3d::Pointer input=reader->GetOutput(); + minc::image3d::Pointer output; + + typedef itk::DanielssonDistanceMapImageFilter< minc::mask3d, minc::image3d > + DistanceMapFilter; + typedef itk::SignedDanielssonDistanceMapImageFilter< minc::mask3d, minc::image3d > + SignedDistanceMapFilter; + + if(ss) + { + SignedDistanceMapFilter::Pointer dist(SignedDistanceMapFilter::New()); + dist->SetInput(input); + dist->Update(); + output=dist->GetOutput(); + } else { + DistanceMapFilter::Pointer dist(DistanceMapFilter::New()); + dist->SetInput(input); + dist->Update(); + output=dist->GetOutput(); + } + + minc::copy_metadata(output,input); + minc::append_history(output,history); + free(history); + + itk::ImageFileWriter< minc::image3d >::Pointer writer = itk::ImageFileWriter<minc::image3d >::New(); + writer->SetFileName(out_f.c_str()); + writer->SetInput( output ); + writer->Update(); + + } catch (const minc::generic_error & err) { + cerr << "Got an error at:" << err.file () << ":" << err.line () << endl; + cerr << err.msg()<<endl; + return 1; + } + return 0; +}
new file mode 100644 --- /dev/null +++ b/minc4itk/tools/itk_resample.cpp @@ -0,0 +1,260 @@ +/* ----------------------------- MNI Header ----------------------------------- +@NAME : itk_resample +@DESCRIPTION: an example of using spline itnerpolation with MINC xfm transform +@COPYRIGHT : + Copyright 2011 Vladimir Fonov, McConnell Brain Imaging Centre, + Montreal Neurological Institute, McGill University. + Permission to use, copy, modify, and distribute this + software and its documentation for any purpose and without + fee is hereby granted, provided that the above copyright + notice appear in all copies. The author and McGill University + make no representations about the suitability of this + software for any purpose. It is provided "as is" without + express or implied warranty. +---------------------------------------------------------------------------- */ +#include <iostream> +#include <fstream> +#include <getopt.h> +#include <vector> +#include <valarray> +#include <math.h> +#include <limits> +#include <unistd.h> + +#include <itkResampleImageFilter.h> +#include <itkAffineTransform.h> +#include <itkNearestNeighborInterpolateImageFunction.h> +#include <itkBSplineInterpolateImageFunction.h> +#include <minc_general_transform.h> + +#include <unistd.h> +#include <getopt.h> +#include <time_stamp.h> // for creating minc style history entry + + +#include <itkImageFileReader.h> +#include <itkImageFileWriter.h> +#include <itkImageIOFactory.h> + +#include "itkMincImageIOFactory.h" +#include "itkMincImageIO.h" +#include "minc_helpers.h" + +//typedef itk::MincImageIO ImageIOType +typedef itk::BSplineInterpolateImageFunction< minc::image3d, double, double > InterpolatorType; +typedef itk::ResampleImageFilter<minc::image3d, minc::image3d> FilterType; +typedef minc::XfmTransform<double,3,3> TransformType; + +using namespace std; + +void show_usage (const char * prog) +{ + std::cerr + << "Usage: "<<prog<<" <input> <output.mnc> " << std::endl + << "--clobber overwrite files" << std::endl + << "--like <example> (default behaviour analogous to use_input_sampling)"<<std::endl + << "--transform <xfm_transform> "<<std::endl + << "--order <n> spline order, default 2 "<<std::endl + << "--uniformize <step> - will make a volume with uniform step size and no direction cosines" << std::endl + << "--invert_transform - apply inverted transform"<<std::endl; +} + +void generate_uniform_sampling(FilterType* flt, const minc::image3d* img,double step) +{ + //obtain physical coordinats of all image corners + minc::image3d::RegionType r=img->GetLargestPossibleRegion(); + std::vector<double> corner[3]; + for(int i=0;i<8;i++) + { + minc::image3d::IndexType idx; + minc::image3d::PointType c; + idx[0]=r.GetIndex()[0]+r.GetSize()[0]*(i%2); + idx[1]=r.GetIndex()[1]+r.GetSize()[1]*((i/2)%2); + idx[2]=r.GetIndex()[2]+r.GetSize()[2]*((i/4)%2); + img->TransformIndexToPhysicalPoint(idx,c); + for(int j=0;j<3;j++) + corner[j].push_back(c[j]); + } + minc::image3d::IndexType start; + FilterType::SizeType size; + FilterType::OriginPointType org; + minc::image3d::SpacingType spc; + spc.Fill(step); + for(int j=0;j<3;j++) + { + std::sort(corner[j].begin(),corner[j].end()); + size[j]=ceil((corner[j][7]-corner[j][0])/step); + org[j]=corner[j][0]; + } + minc::image3d::DirectionType identity; + identity.SetIdentity(); + flt->SetOutputDirection(identity); + start.Fill(0); + flt->SetOutputStartIndex(start); + flt->SetSize(size); + flt->SetOutputOrigin(org); + flt->SetOutputSpacing(spc); +} + +int main (int argc, char **argv) +{ + int verbose=0, clobber=0,skip_grid=0; + int order=2; + std::string like_f,xfm_f,output_f,input_f; + double uniformize=0.0; + int invert=0; + char *history = time_stamp(argc, argv); + + static struct option long_options[] = { + {"verbose", no_argument, &verbose, 1}, + {"quiet", no_argument, &verbose, 0}, + {"clobber", no_argument, &clobber, 1}, + {"like", required_argument, 0, 'l'}, + {"transform", required_argument, 0, 't'}, + {"order", required_argument, 0, 'o'}, + {"uniformize", required_argument, 0, 'u'}, + {"invert_transform", no_argument, &invert, 1}, + {0, 0, 0, 0} + }; + + for (;;) { + /* getopt_long stores the option index here. */ + int option_index = 0; + + int c = getopt_long (argc, argv, "vqcl:t:o:u:", long_options, &option_index); + + /* Detect the end of the options. */ + if (c == -1) break; + + switch (c) + { + case 0: + break; + case 'v': + cout << "Version: 0.1" << endl; + return 0; + case 'l': + like_f=optarg; break; + case 't': + xfm_f=optarg; break; + case 'o': + order=atoi(optarg);break; + case 'u': + uniformize=atof(optarg);break; + case '?': + /* getopt_long already printed an error message. */ + default: + show_usage (argv[0]); + return 1; + } + } + + if ((argc - optind) < 2) { + show_usage(argv[0]); + return 1; + } + input_f=argv[optind]; + output_f=argv[optind+1]; + + if (!clobber && !access (output_f.c_str (), F_OK)) + { + std::cerr << output_f.c_str () << " Exists!" << std::endl; + return 1; + } + + try + { + itk::RegisterMincIO(); + + itk::ImageFileReader<minc::image3d >::Pointer reader = itk::ImageFileReader<minc::image3d >::New(); + + //initializing the reader + reader->SetFileName(input_f.c_str()); + reader->Update(); + + minc::image3d::Pointer in=reader->GetOutput(); + + FilterType::Pointer filter = FilterType::New(); + + //creating coordinate transformation objects + TransformType::Pointer transform = TransformType::New(); + if(!xfm_f.empty()) + { + //reading a minc style xfm file + transform->OpenXfm(xfm_f.c_str()); + if(!invert) transform->Invert(); //should be inverted by default to walk through target space + filter->SetTransform( transform ); + } + + //creating the interpolator + InterpolatorType::Pointer interpolator = InterpolatorType::New(); + interpolator->SetSplineOrder(order); + filter->SetInterpolator( interpolator ); + filter->SetDefaultPixelValue( 0 ); + + //this is for processing using batch system + filter->SetNumberOfThreads(1); + + minc::image3d::Pointer like=0; + if(!like_f.empty()) + { + itk::ImageFileReader<minc::image3d >::Pointer reader = itk::ImageFileReader<minc::image3d >::New(); + reader->SetFileName(like_f.c_str()); + reader->Update(); + if(uniformize!=0.0) + { + generate_uniform_sampling(filter,reader->GetOutput(),uniformize); + } else { + filter->SetOutputParametersFromImage(reader->GetOutput()); + filter->SetOutputDirection(reader->GetOutput()->GetDirection()); + } + like=reader->GetOutput(); + like->DisconnectPipeline(); + } + else + { + if(uniformize!=0.0) + { + generate_uniform_sampling(filter,in,uniformize); + } else { + //we are using original sampling + filter->SetOutputParametersFromImage(in); + filter->SetOutputDirection(in->GetDirection()); + } + } + + filter->SetInput(in); + filter->Update(); + //copy the metadate information, for some reason it is not preserved + //filter->GetOutput()->SetMetaDataDictionary(reader->GetOutput()->GetMetaDataDictionary()); + minc::image3d::Pointer out=filter->GetOutput(); + minc::copy_metadata(out,in); + minc::append_history(out,history); + free(history); + + //correct dimension order + if(like.IsNotNull()) + minc::copy_dimorder(out,like); + + //generic file writer + itk::ImageFileWriter< minc::image3d >::Pointer writer = itk::ImageFileWriter<minc::image3d >::New(); + writer->SetFileName(output_f.c_str()); + + writer->SetInput( out ); + //writer->UseInputMetaDataDictionaryOn(); + + writer->Update(); + + return 0; + } catch (const minc::generic_error & err) { + cerr << "Got an error at:" << err.file () << ":" << err.line () << endl; + return 1; + } + catch( itk::ExceptionObject & err ) + { + std::cerr << "ExceptionObject caught !" << std::endl; + std::cerr << err << std::endl; + return 2; + } + return 0; +};
--- a/progs/CMakeLists.txt +++ b/progs/CMakeLists.txt @@ -155,7 +155,7 @@ DESTINATION bin) # perl and shell scripts -INSTALL(FILES +INSTALL(PROGRAMS mincdiff/mincdiff mincedit/mincedit mincheader/mincheader @@ -163,15 +163,19 @@ PERMISSIONS OWNER_EXECUTE OWNER_WRITE OWNER_READ GROUP_EXECUTE GROUP_READ WORLD_EXECUTE WORLD_READ DESTINATION bin ) +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/minchistory/minchistory.in ${CMAKE_CURRENT_BINARY_DIR}/minchistory @ONLY) +add_custom_target(minchistory chmod +x ${CMAKE_CURRENT_BINARY_DIR}/minchistory DEPENDS ${CMAKE_CURRENT_SOURCE_DIR}/minchistory/minchistory.in) -configure_file(minchistory/minchistory.in ${CMAKE_CURRENT_BINARY_DIR}/minchistory @ONLY) -add_custom_target(minchistory chmod +x ${CMAKE_CURRENT_BINARY_DIR}/minchistory DEPENDS minchistory/minchistory.in) - -configure_file(mincpik/mincpik.in ${CMAKE_CURRENT_BINARY_DIR}/mincpik @ONLY) -add_custom_target(mincpik chmod +x ${CMAKE_CURRENT_BINARY_DIR}/mincpik DEPENDS mincpik/mincpik.in) +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/mincpik/mincpik.in ${CMAKE_CURRENT_BINARY_DIR}/mincpik @ONLY) +add_custom_target(mincpik chmod +x ${CMAKE_CURRENT_BINARY_DIR}/mincpik DEPENDS ${CMAKE_CURRENT_SOURCE_DIR}/mincpik/mincpik.in) INSTALL(FILES ${CMAKE_CURRENT_BINARY_DIR}/minchistory PERMISSIONS OWNER_EXECUTE OWNER_WRITE OWNER_READ GROUP_EXECUTE GROUP_READ WORLD_EXECUTE WORLD_READ DESTINATION bin ) + +INSTALL(PROGRAMS + ${CMAKE_CURRENT_BINARY_DIR}/mincpik + PERMISSIONS OWNER_EXECUTE OWNER_WRITE OWNER_READ GROUP_EXECUTE GROUP_READ WORLD_EXECUTE WORLD_READ + DESTINATION bin )