I accidentally uploaded the segmentation. This is the raw data. The problem occurs with both.

Log output with debug version:

!5.547! WARNING: EnsurePerpendicularNormal(): Plane geometry was _/askew/_, so here it gets rectified by substituting the 3rd column of the indexToWorldMatrix with an appropriate normal vector: [ 0.921861 -0.918454 -0.746066 ], original vector zed was: [ -0.943643 -3.64367e-09 -1.16599 ].
MitkWorkbench: ../../../../MITK/Modules/Core/src/DataManagement/mitkSlicedGeometry3D.cpp:273: virtual void mitk::SlicedGeometry3D::InitializePlanes(const mitk::BaseGeometry*, mitk::PlaneGeometry::PlaneOrientation, bool, bool, bool): Assertion `(standardPlaneNormal - (top ? 1.0 : -1.0) * worldPlaneNormal).GetSquaredNorm() < 0.000001' failed.

@isensee : what is the source of these images? any chance to get the original DICOM?

These images are from the Kora Dataset. Peter Full and @kislinsk have the original DICOM images.
The images were 4d DICOM originally and have been converted into 3d nifti by @kislinsk

@full @schererj @kislinsk : how did you convert these files exactly? If I load the DICOM with a current MITK Workbench -> no problem, re-init fine. Save as nifti, quit, restart, load nifti, re-init -> no problem as well.

So we should investigate two things: 1. what is different in your conversion process? 2. are "your" nifti files somehow invalid and should be rejected on load or should we able to handle them as well.

I converted them with a simple two-liner: mitk::IOUtil::Load<mitk::Image>(); mitk::IOUtil::Save();

Peter and I checked some output NIFTI images and they are valid 3d+t images. Reinit did also work. So for now I suspect that there is another processing step involved? The images above do not match my output files regarding name and file extension I think.

I'll catch some of the files on my way to the Hacking Week...

@full

As far as I can tell files were split along t with fslsplit.

Peter @full converted the files to .nii.gz

I added the _0000 manually. No tool involved here

Just checked the original DICOM and my output file of the derived image above. The original DICOM has the matrix inversion issue already. I bring the files along...

Even when deactivating the askew image handling, there's matrix inversion errors with the original image matrix.

Ideas for further investigation:

• Visualize the axis vectors of the original image matrix (first three columns)
• Check file in other applications (is there a reinit equivalent?)

I have the issue again with other data data from a cooperation partner. Also Slicer has issues with the original DICOM data set:

Geometric issues were found with 1 of the series.  Please use caution.
reg inside examine
Loading with imageIOName: GDCM
Window/level found in DICOM tags (center=212.0, width=487.0) has been applied to volume 25: CINE_segmented_SAX_InlineVF
Irregular volume geometry detected, but maximum error is within tolerance (maximum error of 4.00539e-05 mm, tolerance threshold is 0.001 mm).
Loading with imageIOName: GDCM
Window/level found in DICOM tags (center=212.0, width=487.0) has been applied to volume 25: CINE_segmented_SAX_InlineVF_1
Irregular volume geometry detected, but maximum error is within tolerance (maximum error of 4.00539e-05 mm, tolerance threshold is 0.001 mm).

I guess there's something wrong with the original image matrix.
MITK tries to fix the normal to be perpendicular, which leads to a matrix that can't be inverted anymore.

Looks like the Slicer warning is important, but still it succeeds in loading the data. The numerical error is rather low. Found this discussion which also links to the python code Slicer is using and which emits the warning. Could be of some use for further investigations. https://discourse.slicer.org/t/dicom-scalar-volume-load-irregular-geometry-warning-overly-stringent/3761/14

had a further look into it:

• looked into the matrix vectors. Looks like the first and the second vector is normalized, but the third is not (https://academo.org/demos/3d-vector-plotter/)
• tried to transpose the matrix
• tried to change the last column by sign
• the rectification leads to a normal vector that is the same like the first vector.
• we can set the matrix manually to space directions: (1.9,0,0) (0,1.9,0) (0,0,12) (pixel spacing). But that's not a real solution as we usually can't alter the data.

futher analyses:

• newer data tends to work better
  • newer in terms of date
  • newer in terms of acquisition machine software version
• working and not working data stems from the same machine (Siemens)

In MITK 2016.03, Reinit on the respective data sets work, whereas in MITK 2016.11, it does not work any more.
Hence, it has to do something with the changes introduces by @espak in Geometry (1ca1a7b0574).

could also be that the problems arise due to ITK update (2016.03: 4.7.1, 2016.11: 4.9.0).
We could check that with MITK code from end 2016 and an updated ITK version.

It's apparently an issue of the "dominantAxis" calculation. In some cases this does not make sense and yields the same index twice. Not sure what this "dominantAxis" is supposed to be anyway. Supposedly related to T22254.

The meaning of the dominant axes is also documented here:

https://phabricator.mitk.org/source/mitk/browse/master/Modules/Core/src/DataManagement/mitkPlaneGeometry.cpp$401

It is basically a mapping from the axes of the image coordinate axes to world coordinate axes. It determines which axis of the image coordinate system is the "axial" axis, which one is the "sagittal" and which is the "coronal". Since the "index-to-world" transformation is not just a permutation of the indices and a scaling, but you can rotate the volume in any way, it is not always obvious which axis is which. That's why the "max" call.

If you get the same index twice, that's certainly wrong. It can happen if the plane is rotated by *exactly* 45 degrees around one or two world axes. This case is 'undecided', one choice is not better than the other. You should, however, assign different axes to the different dimensions as 'dominant'. When looking for the dominant axis, you should also check if the value has already been "taken".

Hope that helps.

Does anyone know why the inverse is always calculated and not simply set to the transpose? After all, we are talking about rotation matrices.
In this context, the EnsurePerpendicularNormal method (https://phabricator.mitk.org/source/mitk/browse/master/Modules/Core/src/DataManagement/mitkPlaneGeometry.cpp$37) seems rather complex and also does not catch all cases since it only checks if axis 3 is perpendicular to 1&2 but not if 1&2 are perpendicular to each other. I would prefer to simply check if the matrix R is a rotation matrix by checking

1. det(R) == 1?
2. R*R^T = ID?

I'm also not sure if it is fine to simply substitute one vector if this is not the case. It does not work in all cases and we also do not want to change the data, right? It is further more unclear why this is checked only when calculating the standard planes.

I would propose to check if R is a valid rotation matrix immediately when loading an image and output a warning if this is not the case. For example in cases of minor deviations from condition 2., as it is the case in our test data, it should be up to the user to ignore this. It's a question of how you want to set your epsilon.

Such an "imperfect" rotation matrix in conjunction with "almost" 45° rotations can of course yield the above mentioned error during the dominant axis calculation, but this can easily be caught. I will implement a proposal.

Pushed new branch T25682-DominantAxesError.

While I changed the matrix check, I left it at the current position in in the PlaneGeometry.
We should discuss if it is sensible to move this ckec and also it's applications to a different location, e.g. to SetIndexToWorldTransform of BaseGeometry.

I just realized that it is actually checked in SetIndexToWorldTransform :-D

As I see, the inverse matrix calculation came in with this commit of mine:

https://phabricator.mitk.org/rMITK1b9d9bc730c86c761634ab2a1c5ef042ff32c5c8

which refers to this wiki page from Slicer as 'inspiration':

https://na-mic.org/wiki/Coordinate_System_Conversion_Between_ITK_and_Slicer3

But that page only determines the dominant axes per orientation and does not actually uses inverse matrix calculation.

Most likely I took that code from NifTK that has a viewer that can display several images side-by-side and in image space rather than world space.

https://github.com/NifTK/NifTK/blob/master/MITK/Modules/DnDDisplay/niftkMultiWindowWidget.cxx#L1072

Now I'm curious if it would crash with the attached image.

Another place where inverse matrix calculation has been introduced:

https://phabricator.mitk.org/rMITK113870ce6283059c17436741d05957b276ec22a1