Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

ENH: Added permissive mode to NIFTI reader #4009

Conversation

vfonov
Copy link
Contributor

@vfonov vfonov commented Apr 12, 2023

Added permissive mode to NIFTI reader to allow reading files with non-orthogonal s-form.
Addresses #3994

PR Checklist

  • No API changes were made (or the changes have been approved)
  • No major design changes were made (or the changes have been approved)
  • Added test (or behavior not changed)
  • Updated API documentation (or API not changed)
  • Added license to new files (if any)
  • Added Python wrapping to new files (if any) as described in ITK Software Guide Section 9.5
  • Added ITK examples for all new major features (if any)

Refer to the ITK Software Guide for
further development details if necessary.

@github-actions github-actions bot added type:Infrastructure Infrastructure/ecosystem related changes, such as CMake or buildbots type:Enhancement Improvement of existing methods or implementation type:Testing Ensure that the purpose of a class is met/the results on a wide set of test cases are correct area:IO Issues affecting the IO module type:Data Changes to testing data labels Apr 12, 2023
@dzenanz
Copy link
Member

dzenanz commented Apr 12, 2023

It looks like you need to apply clang-format, see complaints.

CMakeLists.txt Outdated
@@ -338,6 +338,11 @@ if(NOT ("${ITK_NIFTI_IO_ANALYZE_FLAVOR}" STREQUAL "ITK4Warning" OR
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" ON)
mark_as_advanced(ITK_NIFTI_IO_SFORM_PERMISSIVE)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should this be a compile time flag or a runtime environment variable?

Copy link
Member

@issakomi issakomi Apr 12, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It looks like that the option changes the default behavior, not a final behavior:
, m_SFORM_Permissive{ ITK_NIFTI_IO_SFORM_PERMISSIVE }
m_SFORM_Permissive can be changed later itkSetMacro(SFORM_Permissive, bool);

Also I can hardly believe that #cmakedefine ITK_NIFTI_IO_SFORM_PERMISSIVE will work at all.

#cmakedefine VAR ...
will be replaced with either
#define VAR ...
or
/* #undef VAR */

So even if the option is set to "ON", it will not set m_SFORM_Permissive to true here
, m_SFORM_Permissive{ ITK_NIFTI_IO_SFORM_PERMISSIVE },
it will expand to , m_SFORM_Permissive{} and set it to false.
If ITK_NIFTI_IO_SFORM_PERMISSIVE is not defined, the build will fail.

There is #cmakedefine01 s. configure_file, this may work. But, IMHO, it would be better to review the whole thing. The option seems to only set the default, but it is not clear from the description.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

ok, I changed #cmakedefine to #cmakedefine01 and updated option description.

@vfonov vfonov force-pushed the FIX_NIFTI_NON_ORTHO_COORDINATES branch 2 times, most recently from c745338 to 0b9b6e0 Compare April 12, 2023 23:00
CMakeLists.txt Outdated
@@ -338,6 +338,11 @@ if(NOT ("${ITK_NIFTI_IO_ANALYZE_FLAVOR}" STREQUAL "ITK4Warning" OR
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)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am very opposed to this being the default compile time option.

Was OK with a runtime selection by the developer to explicitly request to ignore the fact that this behavior is not conformant with the rest of the ITK standard.

I think it is a very dangerous proposition to allow very hard to understand data assumption violations as the default.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why is a compile time option needed? The runtime behavior should be backward compatible with the historical behavior, and developers who wish to change the behavior have the option to set the desired new behavior at runtime.

PS: Thanks for putting this patch forward for discussion.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So, I thought this way we will restore compatibility with the version before the strict checking was introduced.

Copy link
Member

@blowekamp blowekamp Apr 13, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am also apposed to a compile time option.

Perhaps a "GlobalDefault" runtime variable would be more appropriate? With the default as Hans has described.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

so, what's the current recommendation?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

ITK_NIFTI_IO_SFORM_PERMISSIVE defaulting to OFF? And whoever doesn't like that default, can change it at both configure time and run time.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Added permissive mode to NIFTI reader to allow reading
files with non-orthogonal s-form.
@vfonov vfonov force-pushed the FIX_NIFTI_NON_ORTHO_COORDINATES branch from 0b9b6e0 to 3d3f3e5 Compare April 13, 2023 14:52
Made ITK_NIFTI_IO_SFORM_PERMISSIVE_DEFAULT OFF by default,
added support for environment variable ITK_NIFTI_SFORM_PERMISSIVE
@vfonov vfonov force-pushed the FIX_NIFTI_NON_ORTHO_COORDINATES branch from be4c358 to 3e45ceb Compare April 13, 2023 16:57
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]));
Copy link
Member

@issakomi issakomi Apr 14, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

itksys::SystemTools::PutEnv("ITK_NIFTI_SFORM_PERMISSIVE"); doesn't work

~$ export ITK_NIFTI_SFORM_PERMISSIVE
~$ env | grep ITK_NIFTI_SFORM_PERMISSIVE
~$

Probably at least empty value is required for environment variable.
Test doesn't fail without 2 above lines (174 and 175).

BTW, FYI, another test itkNiftiReadWriteDirectionSmallVoxelTest fails because number of arguments is wrong now.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Export and PutEnv might work differently?

Copy link
Member

@issakomi issakomi Apr 14, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Export and PutEnv might work differently?

Yes, export doesn't override last value if there is no "="

r@deb2:~$ export SOME_VAR=1
r@deb2:~$ export SOME_VAR
r@deb2:~$ env | grep SOME_VAR
SOME_VAR=1
r@deb2:~$ 

PutEnv seems to unset the variable if there is no "=". Anyway, the test works without these 2 lines (174, 175).

bool SystemTools::PutEnv(const std::string& env)
{
  size_t pos = env.find('=');
  if (pos != std::string::npos) {
    std::string name = env.substr(0, pos);
    return setenv(name.c_str(), env.c_str() + pos + 1, 1) == 0;
  } else {
    return kwsysUnPutEnv(env) == 0;
  }
}

// Fix to deal with non-orthogonal matrixes
// this approach uses SVD decomposition
// maybe it is better to use nifti_make_orthog_mat44
// to be consistent with other software?
Copy link
Member

@issakomi issakomi Apr 16, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Unfortunately the block below this line doesn't produce orthogonal matrix. The test image NonOrthoSform.nii.gz ends up with ITK image:

  Direction: 
0.694449 0.264737 0.469327
-0.625039 -0.670148 -0.125412
0.356465 -0.693409 0.874073

The matrix is not orthogonal. I have tested the alternative you mentioned above and it works very good, nifti_make_orthog_mat44 generates orthogonal matrices.

Updated post!

Also the correct origin is [-120.117, 77.8649, -128.697] , not [-0.724625, -0.614831, 0.311127].

Probably this works:

      else if (m_SFORM_Permissive)
      {
        // nifti_make_orthog_mat44 does 3x3 orthogonal matrix, other components are 0 and [3][3] is 1
        mat44 _mat = nifti_make_orthog_mat44(this->m_NiftiImage->sto_xyz.m[0][0],
                                             this->m_NiftiImage->sto_xyz.m[0][1],
                                             this->m_NiftiImage->sto_xyz.m[0][2],
                                             this->m_NiftiImage->sto_xyz.m[1][0],
                                             this->m_NiftiImage->sto_xyz.m[1][1],
                                             this->m_NiftiImage->sto_xyz.m[1][2],
                                             this->m_NiftiImage->sto_xyz.m[2][0],
                                             this->m_NiftiImage->sto_xyz.m[2][1],
                                             this->m_NiftiImage->sto_xyz.m[2][2]);
        _mat.m[0][3] = this->m_NiftiImage->sto_xyz.m[0][3]; // Origin X
        _mat.m[1][3] = this->m_NiftiImage->sto_xyz.m[1][3]; // Origin Y
        _mat.m[2][3] = this->m_NiftiImage->sto_xyz.m[2][3]; // Origin Z
        // Probably not required to copy?
        /*
        _mat.m[3][3] = this->m_NiftiImage->sto_xyz.m[3][3]; // 1
        _mat.m[3][0] = this->m_NiftiImage->sto_xyz.m[3][0]; // 0
        _mat.m[3][1] = this->m_NiftiImage->sto_xyz.m[3][1]; // 0
        _mat.m[3][2] = this->m_NiftiImage->sto_xyz.m[3][2]; // 0
        */
        // Copying scales seems to be not required for spacing, it is not taken from the matrix
        /*
        const vnl_matrix_fixed<float, 4, 4> sform_as_matrix{ &(this->m_NiftiImage->sto_xyz.m[0][0]) };
        const vnl_matrix_fixed<float, 3, 3> rotation = sform_as_matrix.extract(3, 3, 0, 0);
        _mat.m[0][0] *= rotation.get_column(0).magnitude();
        _mat.m[1][0] *= rotation.get_column(0).magnitude();
        _mat.m[2][0] *= rotation.get_column(0).magnitude();
        _mat.m[0][1] *= rotation.get_column(1).magnitude();
        _mat.m[1][1] *= rotation.get_column(1).magnitude();
        _mat.m[2][1] *= rotation.get_column(1).magnitude();
        _mat.m[0][1] *= rotation.get_column(2).magnitude();
        _mat.m[1][1] *= rotation.get_column(2).magnitude();
        _mat.m[2][1] *= rotation.get_column(2).magnitude();
        */
        itkWarningMacro(<< this->GetFileName() << " has non-orthogonal sform, corrected");
        this->m_SFORM_Corrected = true;
        return _mat;
      }

It is still in testing!

The current code takes the whole 4x4 matrix, but 4th column is an origin (translation)
mat = svd.V() * svd.U().conjugate_transpose();, it is probably the issue.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It appears that the logic for making an orthogonal matrix is not correct. See #4009 (comment)

OK, it was my mistake - I should have done it only on 3x3 upper left submatrix. Fixed now.

Copy link
Member

@hjmjohnson hjmjohnson left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It appears that the logic for making an orthogonal matrix is not correct. See #4009 (comment)

@vfonov vfonov requested a review from hjmjohnson April 17, 2023 22:53
@hjmjohnson
Copy link
Member

@issakomi Would you please make a final review of the latest changes to ensure that they address your previously identified bug?

Thanks,
Hans

@issakomi
Copy link
Member

@issakomi Would you please make a final review of the latest changes to ensure that they address your previously identified bug?

The bug is fixed. Thanks.

} // sform is orthonormal
} // sform not NIFTI_XFORM_UNKNOWN
}
else if (m_SFORM_Permissive) // sform is orthonormal
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Very minor comment. IMHO, the comment // sform is orthonormal marks the end of the previous block. At this line it might be a little bit misleading, IMHO.

@dzenanz dzenanz merged commit 815d0f0 into InsightSoftwareConsortium:master Apr 20, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
area:IO Issues affecting the IO module type:Data Changes to testing data type:Enhancement Improvement of existing methods or implementation type:Infrastructure Infrastructure/ecosystem related changes, such as CMake or buildbots type:Testing Ensure that the purpose of a class is met/the results on a wide set of test cases are correct
Projects
None yet
Development

Successfully merging this pull request may close these issues.

5 participants