Skip to content

Commit

Permalink
ENH: Made ITK_NIFTI_IO_SFORM_PERMISSIVE_DEFAULT OFF by default
Browse files Browse the repository at this point in the history
Made ITK_NIFTI_IO_SFORM_PERMISSIVE_DEFAULT OFF by default,
added support for environment variable ITK_NIFTI_SFORM_PERMISSIVE
  • Loading branch information
vfonov committed Apr 13, 2023
1 parent 3d3f3e5 commit 3e45ceb
Show file tree
Hide file tree
Showing 4 changed files with 25 additions and 4 deletions.
4 changes: 2 additions & 2 deletions CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -339,8 +339,8 @@ endif()
mark_as_advanced(ITK_NIFTI_IO_ANALYZE_FLAVOR)

# Allow non-orthogonal rotation matrix approximation in NIFTI files
option(ITK_NIFTI_IO_SFORM_PERMISSIVE "Allow non-orthogonal rotation matrix in NIFTI sform by default" ON)
mark_as_advanced(ITK_NIFTI_IO_SFORM_PERMISSIVE)
option(ITK_NIFTI_IO_SFORM_PERMISSIVE_DEFAULT "Allow non-orthogonal rotation matrix in NIFTI sform by default" OFF)
mark_as_advanced(ITK_NIFTI_IO_SFORM_PERMISSIVE_DEFAULT)



Expand Down
10 changes: 9 additions & 1 deletion Modules/IO/NIFTI/src/itkNiftiImageIO.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -23,6 +23,8 @@
#include "itkNiftiImageIOConfigurePrivate.h"
#include "itkMakeUniqueForOverwrite.h"
#include "vnl/algo/vnl_svd.h"
#include "itksys/SystemTools.hxx"
#include "itksys/SystemInformation.hxx"

namespace itk
{
Expand Down Expand Up @@ -450,7 +452,7 @@ NiftiImageIO::NiftiImageIO()
: m_NiftiImageHolder(new NiftiImageProxy(nullptr))
, m_NiftiImage(*m_NiftiImageHolder.get())
, m_LegacyAnalyze75Mode{ ITK_NIFTI_IO_ANALYZE_FLAVOR_DEFAULT }
, m_SFORM_Permissive{ ITK_NIFTI_IO_SFORM_PERMISSIVE }
, m_SFORM_Permissive{ ITK_NIFTI_IO_SFORM_PERMISSIVE_DEFAULT }
{
this->SetNumberOfDimensions(3);
nifti_set_debug_level(0); // suppress error messages
Expand All @@ -462,6 +464,12 @@ NiftiImageIO::NiftiImageIO()
this->AddSupportedWriteExtension(ext);
this->AddSupportedReadExtension(ext);
}
std::string envVar;
if (itksys::SystemTools::GetEnv("ITK_NIFTI_SFORM_PERMISSIVE", envVar))
{
envVar = itksys::SystemTools::UpperCase(envVar);
this->SetSFORM_Permissive(envVar != "NO" && envVar != "OFF" && envVar != "FALSE");
}
}

NiftiImageIO::~NiftiImageIO()
Expand Down
2 changes: 1 addition & 1 deletion Modules/IO/NIFTI/src/itkNiftiImageIOConfigurePrivate.h.in
Original file line number Diff line number Diff line change
Expand Up @@ -24,6 +24,6 @@
// Analyze default flavour
#define ITK_NIFTI_IO_ANALYZE_FLAVOR_DEFAULT NiftiImageIOEnums::Analyze75Flavor::Analyze@ITK_NIFTI_IO_ANALYZE_FLAVOR@
// Enable permissive treatment of SFORM
#cmakedefine01 ITK_NIFTI_IO_SFORM_PERMISSIVE
#cmakedefine01 ITK_NIFTI_IO_SFORM_PERMISSIVE_DEFAULT

#endif // itkNiftiImageIOConfigurePrivate_h
13 changes: 13 additions & 0 deletions Modules/IO/NIFTI/test/itkNiftiReadWriteDirectionTest.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,8 @@
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"
#include "itkTestingMacros.h"
#include "itksys/SystemTools.hxx"
#include "itksys/SystemInformation.hxx"

// debug
#include <map>
Expand Down Expand Up @@ -161,6 +163,17 @@ itkNiftiReadWriteDirectionTest(int argc, char * argv[])
// should not read the image with non-orthogonal sform
ITK_TRY_EXPECT_EXCEPTION(ReadImage<TestImageType>(argv[4], false));

// check if environment flag is used properly
// check without permissive option
itksys::SystemTools::PutEnv("ITK_NIFTI_SFORM_PERMISSIVE=NO");
ITK_TRY_EXPECT_EXCEPTION(itk::ReadImage<TestImageType>(argv[4]));
// check with permissive option
itksys::SystemTools::PutEnv("ITK_NIFTI_SFORM_PERMISSIVE=YES");
ITK_TRY_EXPECT_NO_EXCEPTION(itk::ReadImage<TestImageType>(argv[4]));
// check with permissive option, when nothing is specified
itksys::SystemTools::PutEnv("ITK_NIFTI_SFORM_PERMISSIVE");
ITK_TRY_EXPECT_NO_EXCEPTION(itk::ReadImage<TestImageType>(argv[4]));

// this should work
TestImageType::Pointer inputImageNonOrthoSform = ReadImage<TestImageType>(argv[4], true);
dictionary = inputImageNonOrthoSform->GetMetaDataDictionary();
Expand Down

0 comments on commit 3e45ceb

Please sign in to comment.