FEATURE 609 Save of volumic images. Anatomical orientation, translation...

FEATURE 609 Save of volumic images. Anatomical orientation, translation (Offset) and rotation (TransformMatrix) taken into account. 
Still some works need to be done : 
* Action move frame should not work on frame coming from reorientation => Improve action ergonomy, mainly + define special frame for reorientation.
* Saving images in MHA has some weird behavior since we open, modify .mha, and this one contains raw image data => more test or decide to avoid .mha file format and prefer using .mhd (better).
* Dicom component does not (at the moment) read rotation information from tag => Improve it to read the data from the tag
* Saving image NOT in RAI, after user modification of the frame becomes WRONG => decide how to handle 3D transformation and frame reorientation. We should separate them.
* After saving (save as), input image data changes, its frame is reset => Need to reload correctly the image (ImageComponent::postSaving() method).

git-svn-id: svn+ssh://scm.forge.imag.fr/var/lib/gforge/chroot/scmrepos/svn/camitk/trunk/camitk@1963 ec899d31-69d1-42ba-9299-647d76f65fb3

imaging/components/dicom/DicomComponentExtension.cpp

@@ -44,7 +44,7 @@
 using namespace camitk;
 
 // --------------- Extension Declaration -------------------
-Q_EXPORT_PLUGIN2(dicomcomponentextension, DicomComponentExtension);
+Q_EXPORT_PLUGIN2(dicomcomponentextension, DicomComponentExtension)
 
 // --------------- GetFileExtensions -------------------
 QStringList DicomComponentExtension::getFileExtensions() const {

imaging/components/dicom/DicomComponentExtension.h

@@ -30,8 +30,9 @@
 
 #include<ComponentExtension.h>
 
-namespace camitk {
-class Component;
+namespace camitk
+{
+    class Component;
 }
 
 class DicomDialog;

sdk/actions/image/reorientimage/ReorientImage.cpp

@@ -328,7 +328,7 @@ Action::ApplyStatus ReorientImage::process( ImageComponent * image) {
         vtkSmartPointer<vtkMatrix4x4> originalMatrix = image->getTransform()->GetMatrix();
 
         if (ui.resetImageOriginRadioButton->isChecked()) {
-            vtkSmartPointer<vtkMatrix4x4> backToOriginMatrix = ImageOrientationHelper::getTransformFromRAI(orient, dims[0], dims[1], dims[2]);
+            vtkSmartPointer<vtkMatrix4x4> backToOriginMatrix = ImageOrientationHelper::getTransformToRAI(orient, dims[0], dims[1], dims[2]);
             vtkSmartPointer<vtkTransform> backToOrigin = vtkSmartPointer<vtkTransform>::New();
             backToOrigin->SetMatrix(backToOriginMatrix);
             image->setTransform(backToOrigin);
@@ -355,7 +355,7 @@ Action::ApplyStatus ReorientImage::process( ImageComponent * image) {
 
         // Move the image back to the origin thanks to Frame
         if (ui.resetImageOriginRadioButton->isChecked()) {
-            vtkSmartPointer<vtkMatrix4x4> backToOriginMatrix = ImageOrientationHelper::getTransformFromRAI(orient, dims[0], dims[1], dims[2]);
+            vtkSmartPointer<vtkMatrix4x4> backToOriginMatrix = ImageOrientationHelper::getTransformToRAI(orient, dims[0], dims[1], dims[2]);
             vtkSmartPointer<vtkTransform> backToOrigin = vtkSmartPointer<vtkTransform>::New();
             backToOrigin->SetMatrix(backToOriginMatrix);
             result->setTransform(backToOrigin);

sdk/components/vtkimage/VtkImageComponent.cpp

@@ -28,6 +28,8 @@
 
 //-- Qt
 #include <QFileInfo>
+#include <QRegExp>
+#include <QTextStream>
 
 //-- Vtk
 #include <vtkImageData.h>
@@ -40,6 +42,7 @@
 #include <vtkPNMReader.h>
 #include <vtkMetaImageReader.h>
 #include <vtkImageFlip.h>
+#include <vtkInformation.h>
 
 using namespace camitk;
 
@@ -95,12 +98,15 @@ void VtkImageComponent::createComponent(const QString& filename) {
                 image = reader->GetOutput();
             }
 
-            // Get the image orientation when possible
+            // Get the image orientation & rotation when possible
             ImageOrientationHelper::PossibleImageOrientations orientation = ImageOrientationHelper::RAI;
-            if(vtkSmartPointer<vtkMetaImageReader> metaImageReader = vtkMetaImageReader::SafeDownCast(reader))
+            vtkSmartPointer<vtkMatrix4x4> rotationMatrix;
+            if(vtkSmartPointer<vtkMetaImageReader> metaImageReader = vtkMetaImageReader::SafeDownCast(reader)){
                 orientation = ImageOrientationHelper::getOrientationAsEnum(QString(metaImageReader->GetAnatomicalOrientation()));
+                rotationMatrix = readMetaImageTransformMatrix(filename);
+            }
 
-            setImageData(image, false, orientation);
+            setImageData(image, false, orientation, rotationMatrix);
         }
 
         catch (AbortException e) {
@@ -108,3 +114,36 @@ void VtkImageComponent::createComponent(const QString& filename) {
         }
     }
 }
+
+// -------------------- readMetaImageTransformMatrix --------------------
+vtkSmartPointer<vtkMatrix4x4> VtkImageComponent::readMetaImageTransformMatrix(const QString & fileName){
+
+    QFile file(fileName);
+    if (!file.open(QIODevice::ReadOnly | QIODevice::Text))
+        CAMITK_ERROR("VtkImageComponent", "readMetaImageTransformMatrix", "Error reading TransformMatrix tag from image " + fileName.toStdString())
+
+    QTextStream in(&file);
+    QString line = in.readLine();
+    while (!line.isNull()) {
+        // Find lind feature TransformMatrix tag information
+        QRegExp regExp = QRegExp("TransformMatrix = ([-+]?[0-9]*\\.?[0-9]+([eE][-+]?[0-9]+)?[ |\\n])+");
+        if(line.contains(regExp)){ // found it !
+            // retrieve information
+            QStringList values = line.split(" ");
+            vtkSmartPointer<vtkMatrix4x4> rotationMatrix = vtkMatrix4x4::New();
+            for(int j=0; j < 3; j++){
+                for(int i=0; i < 3; i++){
+                    rotationMatrix->SetElement(i, j, values.at(j*3 + i + 2).toDouble());
+                }
+                std::cout << std::endl;
+            }
+            file.close();
+
+            return rotationMatrix;
+        }
+        line = in.readLine();
+    }
+
+    file.close(); // close the file handle.
+    return NULL;
+}

sdk/components/vtkimage/VtkImageComponent.h

@@ -46,10 +46,18 @@ public:
     VtkImageComponent(const QString&) throw (camitk::AbortException);
 
     /// needed for deleting
-    virtual ~VtkImageComponent() {};
+    virtual ~VtkImageComponent() {}
 
     /// actually create the component from the file
     virtual void createComponent(const QString&);
+
+private:
+    /**
+     * Read the TranformMatrix tag (rotation of the image) from the input MHA/MHD image
+     * @param fileName The input MHA / MHD file name to read the tag from.
+     * @return A 4x4 homogenious matrix containing the Rotation information of the image
+     */
+    vtkSmartPointer<vtkMatrix4x4> readMetaImageTransformMatrix(const QString & fileName);
 };
 
 

sdk/components/vtkimage/VtkImageComponentExtension.cpp

@@ -45,7 +45,7 @@ using namespace camitk;
 #include <vtkMetaImageWriter.h>
 
 // --------------- declare the plugin -------------------
-Q_EXPORT_PLUGIN2(vtkimage, VtkImageComponentExtension);
+Q_EXPORT_PLUGIN2(vtkimage, VtkImageComponentExtension)
 
 
 // --------------- getName -------------------
@@ -157,7 +157,9 @@ bool VtkImageComponentExtension::save(Component* component) const {
                 writer->SetFileDimensionality(2);
             } else if ((QString::compare(fileExt, "mhd", Qt::CaseInsensitive) == 0) ||
                        (QString::compare(fileExt, "mha", Qt::CaseInsensitive) == 0)) {
-                writer = vtkSmartPointer<vtkMetaImageWriter>::New();
+                vtkSmartPointer<vtkMetaImageWriter> metaImgWriter = vtkSmartPointer<vtkMetaImageWriter>::New();
+                metaImgWriter->SetCompression(false);
+                writer = metaImgWriter;
                 writer->SetFileDimensionality(3);
             } else {
                 QMessageBox::warning(NULL, "Saving Error", tr("Cannot save file: unrecognized extension \".") + fileExt + "\"");
@@ -179,7 +181,7 @@ bool VtkImageComponentExtension::save(Component* component) const {
                 // vtkImageData or vtkMetaImageWriter classes
                 // we are obliged to open the whole image, as a text, replace the anatomical orientation string
                 // store it back in ASCII (to preserve data format)
-                if (QString::compare(fileExt, "mha", Qt::CaseInsensitive) == 0) {
+                if ((QString::compare(fileExt, "mha", Qt::CaseInsensitive) == 0) || (QString::compare(fileExt, "mhd", Qt::CaseInsensitive) == 0)) {
                     QByteArray fileData;
                     QFile file(component->getFileName());
                     file.open(QIODevice::ReadWrite);
@@ -189,35 +191,22 @@ bool VtkImageComponentExtension::save(Component* component) const {
                     // Store anatomical orientation in header
                     text.replace(QString("AnatomicalOrientation = ???"), QString("AnatomicalOrientation = ") + ImageOrientationHelper::getOrientationAsQString(img->getInitialOrientation()));
 
-                    file.seek(0); // go to the beginning of the file
-                    file.write(text.toAscii()); // write back the file, in ASCII to preserve data format
-
-                    file.close(); // close the file handle.
-                }
-
-                if (QString::compare(fileExt, "mhd", Qt::CaseInsensitive) == 0) {
-                    QByteArray fileData;
-                    QFile file(component->getFileName());
-                    file.open(QIODevice::ReadWrite);
-                    fileData = file.readAll();
-                    // empty the file
-                    file.remove();
-                    file.open(QIODevice::ReadWrite);
-
-                    QString text(fileData);
-                    // Store anatomical orientation in header
-                    text.replace(QString("AnatomicalOrientation = ???"), QString("AnatomicalOrientation = ") + ImageOrientationHelper::getOrientationAsQString(img->getInitialOrientation()));
-                    // Remove transform matrix, offset and center of orientation tag (added by VTK and useless for us).
-//                    text.replace(QRegExp("TransformMatrix = (([-+]?(?:\\b[0-9]+(?:\\.[0-9]*)?|\\.[0-9]+\\b)(?:[eE][-+]?[0-9]+\\b)?)[ |\\n])+"), QString(""));
-//                    text.replace(QRegExp("Offset = (([-+]?(?:\\b[0-9]+(?:\\.[0-9]*)?|\\.[0-9]+\\b)(?:[eE][-+]?[0-9]+\\b)?)[ |\\n])+"), QString(""));
-//                    text.replace(QRegExp("CenterOfRotation = (([-+]?(?:\\b[0-9]+(?:\\.[0-9]*)?|\\.[0-9]+\\b)(?:[eE][-+]?[0-9]+\\b)?)[ |\\n])+"), QString(""));
+                    // Retrieve orientation information
+                    vtkSmartPointer<vtkMatrix4x4> rotationMatrix = img->getRotationMatrix();
+                    QString rotationInfo = "TransformMatrix = ";
+                    for(int j=0; j < 3; j++){
+                        for(int i=0; i < 3; i++){
+                            rotationInfo += QString::number(rotationMatrix->GetElement(i, j)) + " ";
+                        }
+                    }
+                    rotationInfo += "\n";
+                    text.replace(QRegExp("TransformMatrix = (([-+]?(?:\\b[0-9]+(?:\\.[0-9]*)?|\\.[0-9]+\\b)(?:[eE][-+]?[0-9]+\\b)?)[ |\\n])+"), rotationInfo);
 
                     file.seek(0); // go to the beginning of the file
                     file.write(text.toAscii()); // write back the file, in ASCII to preserve data format
 
                     file.close(); // close the file handle.
                 }
-
             }
 
             return true;

sdk/libraries/core/component/image/ImageComponent.cpp

@@ -219,9 +219,6 @@ void ImageComponent::setImageData(vtkSmartPointer<vtkImageData> anImageData,
                                   ImageOrientationHelper::PossibleImageOrientations initialOrientation,
                                   vtkSmartPointer<vtkMatrix4x4> initialRotationMatrix ) {
 
-    double* orig = anImageData->GetOrigin();
-    CAMITK_DEBUG("ImageComponent", "setImageData", "initial origin = (" + QString::number(orig[0]).toStdString() + ", " + QString::number(orig[1]).toStdString() + ", " + QString::number(orig[2]).toStdString() + ")")
-
     this->initialOrientation = initialOrientation;
 
     originalImageData = NULL;
@@ -239,21 +236,13 @@ void ImageComponent::setImageData(vtkSmartPointer<vtkImageData> anImageData,
     } else
         originalImageData = anImageData;
 
-    // 1. Get / compute the initial transformation of the image
-    vtkSmartPointer<vtkMatrix4x4> initialTransformationMatrix = vtkSmartPointer<vtkMatrix4x4>::New();
-    // a. Get translation from the origin of the vtkImageData
+    // 1. Get / compute the initial translation of the image
     double t_x, t_y, t_z;
     originalImageData->GetOrigin(t_x, t_y, t_z);
     vtkSmartPointer<vtkTransform> initialTranslation = vtkTransform::New();
     initialTranslation->Identity();
     initialTranslation->Translate(t_x, t_y, t_z);
     initialTranslation->Update();
-    // b. Get the rotation if provided
-    if(initialRotationMatrix){
-        vtkMatrix4x4::Multiply4x4(initialTranslation->GetMatrix(), initialRotationMatrix, initialTransformationMatrix);
-    }else{
-        initialTransformationMatrix = initialTranslation->GetMatrix();
-    }
 
     // 2. Get the orientation of the image
     double* imgSpacing = originalImageData->GetSpacing();
@@ -269,33 +258,37 @@ void ImageComponent::setImageData(vtkSmartPointer<vtkImageData> anImageData,
     orientationProperty->setReadOnly(true);
     addProperty(orientationProperty);
 
-    // 4. Compute the full transform (translation + rotation + orientation)
-    vtkSmartPointer<vtkMatrix4x4> transformMatrix = vtkSmartPointer<vtkMatrix4x4>::New();
-    vtkMatrix4x4::Multiply4x4(initialTransformationMatrix, orientationToRAIMatrix, transformMatrix);
+    // 4. Store the Translation * orientation -> RAI = initialImageDataTransform
+    initialImageDataTransform = vtkTransform::New();
+    vtkSmartPointer<vtkMatrix4x4> initialImageDataMatrix = vtkMatrix4x4::New();
+    vtkMatrix4x4::Multiply4x4(initialTranslation->GetMatrix(), orientationToRAIMatrix, initialImageDataMatrix);
+    initialImageDataTransform->SetMatrix(initialImageDataMatrix);
+
+    // 5. Store the rotation * translation * orientation -> RAI = initialFrameTransform
+    initialFrameTransform = vtkTransform::New();
+    vtkSmartPointer<vtkMatrix4x4> initialFrameMatrix = vtkMatrix4x4::New();
+    if(!initialRotationMatrix){
+        initialRotationMatrix = vtkMatrix4x4::New();
+        initialRotationMatrix->Identity();
+    }
+    rotationMatrix = initialRotationMatrix;
+    vtkMatrix4x4::Multiply4x4(initialRotationMatrix, initialTranslation->GetMatrix(), initialFrameMatrix);
+    vtkMatrix4x4::Multiply4x4(initialFrameMatrix, orientationToRAIMatrix, initialFrameMatrix);
+    initialFrameTransform->SetMatrix(initialFrameMatrix);
 
-    // 5. Apply this transform to the data (when using rotation in a Matrix -> need to use vtkImageReslice filter)
+    // 6. Apply initialImageDataTransform to the image data
     vtkSmartPointer<vtkImageReslice> imageResliceFilter = vtkSmartPointer<vtkImageReslice>::New();
     imageResliceFilter->SetInput(anImageData);
     imageResliceFilter->SetOutputDimensionality(3);
-    // imageResliceFilter->SetOutputOrigin(0, 0, 0);
-    imageResliceFilter->SetResliceAxes(transformMatrix);
+    imageResliceFilter->SetResliceAxes(initialImageDataTransform->GetMatrix());
     originalImageData = imageResliceFilter->GetOutput();
     originalImageData->Update();
 
-    // 6. Compute the inverse transformation
-    vtkSmartPointer<vtkTransform> inverseTransform = vtkTransform::New();
-    inverseTransform->SetMatrix(transformMatrix);
-    inverseTransform->Inverse();
-    inverseTransform->Update();
-
-    // 7. Save this inverse transform in case of saving
-    backToOriginalTransform = inverseTransform;
-
-    // 8. Store this transform as the current frame
+    // 7. Store initialFrameTransform as the current frame
     // note: we need to get another matrix instance for the transformation (else it is deleted)
     vtkSmartPointer<vtkMatrix4x4> originalMatrix = this->getTransform()->GetMatrix();
     vtkSmartPointer<vtkMatrix4x4> updatedMatrix = vtkSmartPointer<vtkMatrix4x4>::New();
-    vtkMatrix4x4::Multiply4x4(originalMatrix, transformMatrix, updatedMatrix);
+    vtkMatrix4x4::Multiply4x4(originalMatrix, initialFrameTransform->GetMatrix(), updatedMatrix);
     vtkSmartPointer<vtkTransform> updatedTransform = vtkSmartPointer<vtkTransform>::New();
     updatedTransform->SetMatrix(updatedMatrix);
     setTransform(updatedTransform);
@@ -308,37 +301,61 @@ void ImageComponent::setImageData(vtkSmartPointer<vtkImageData> anImageData,
 
     // Init or update camitk:Properties
     updateImageProperties();
-
-    orig = originalImageData->GetOrigin();
-    CAMITK_DEBUG("ImageComponent", "setImageData", "end initial origin = (" + QString::number(orig[0]).toStdString() + ", " + QString::number(orig[1]).toStdString() + ", " + QString::number(orig[2]).toStdString() + ")")
-
 }
 
 // -------------------- prepareForSaving --------------------
 void ImageComponent::prepareForSaving() {
 
-    // 1st apply inverse orientation transform
-    vtkSmartPointer<vtkImageReslice> backToOrigingFilter = vtkSmartPointer<vtkImageReslice>::New();
-    backToOrigingFilter->SetInput(originalImageData);
-    backToOrigingFilter->SetOutputDimensionality(3);
-    backToOrigingFilter->SetResliceAxes(backToOriginalTransform->GetMatrix());
-    originalImageData = backToOrigingFilter->GetOutput();
+    // 1. Put back image data to its original state
+    vtkSmartPointer<vtkMatrix4x4> backToOriginalImageDataMatrix = vtkMatrix4x4::New();
+    backToOriginalImageDataMatrix->DeepCopy(initialImageDataTransform->GetMatrix());
+    backToOriginalImageDataMatrix->Invert();
+    vtkSmartPointer<vtkImageReslice> backToOrigingImageDataFilter = vtkSmartPointer<vtkImageReslice>::New();
+    backToOrigingImageDataFilter->SetInput(originalImageData);
+    backToOrigingImageDataFilter->SetOutputDimensionality(3);
+    backToOrigingImageDataFilter->SetResliceAxes(backToOriginalImageDataMatrix);
+    originalImageData = backToOrigingImageDataFilter->GetOutput();
     originalImageData->Update();
 
-    // 2. Retrieve the user transformation made since opening of the image by cancelling
-    // the transformations made by opening the image stored in the frame.
+    // 2. Retrieve the matrice of the user Muser which has been applied to the initial frame
+    // Assuming Mframe = MinitFrame * Muser
+    vtkSmartPointer<vtkMatrix4x4> initialFrameMatrixInverse = vtkMatrix4x4::New();
+    initialFrameMatrixInverse->DeepCopy(initialFrameTransform->GetMatrix());
+    initialFrameMatrixInverse->Invert();
     const vtkSmartPointer<vtkTransform> frame = getTransform();
     vtkSmartPointer<vtkMatrix4x4> frameMatrix = frame->GetMatrix();
     vtkSmartPointer<vtkMatrix4x4> userMatrix = vtkMatrix4x4::New();
-    vtkMatrix4x4::Multiply4x4(backToOriginalTransform->GetMatrix(), frameMatrix, userMatrix);
-
-    // 3. Apply user modifications to the image data to save them.
-    vtkSmartPointer<vtkImageReslice> userFilter = vtkSmartPointer<vtkImageReslice>::New();
-    userFilter->SetInput(originalImageData);
-    userFilter->SetOutputDimensionality(3);
-    userFilter->SetResliceAxes(userMatrix);
-    originalImageData = userFilter->GetOutput();
+    vtkMatrix4x4::Multiply4x4(initialFrameMatrixInverse, frameMatrix , userMatrix);
+    vtkSmartPointer<vtkTransform> userTransform = vtkTransform::New();
+    userTransform->SetMatrix(userMatrix);
+
+    // 3. Retrieve the translation of the user matrix and apply it to the image data
+    // This way, the translation will be saved under the Offset tag
+    double pos[3];
+    userTransform->GetPosition(pos);
+    vtkSmartPointer<vtkTransform> translationTransform = vtkTransform::New();
+    translationTransform->Identity();
+    translationTransform->Translate(pos[0], pos[1], pos[2]);
+    translationTransform->Update();
+    vtkSmartPointer<vtkImageReslice> translationFilter = vtkSmartPointer<vtkImageReslice>::New();
+    translationFilter->SetInput(originalImageData);
+    translationFilter->SetOutputDimensionality(3);
+    translationFilter->SetResliceAxes(translationTransform->GetMatrix());
+    originalImageData = translationFilter->GetOutput();
     originalImageData->Update();
+
+    // 4. Retrieve the rotation matrix from the user matrix to manually save it
+    // only if user has modify it
+    if( (userMatrix->GetElement(0,0) != 1) ||
+        (userMatrix->GetElement(1,1) != 1) ||
+        (userMatrix->GetElement(2,2) != 1) ){
+
+        rotationMatrix = vtkMatrix4x4::New();
+        rotationMatrix->DeepCopy(userMatrix);
+        rotationMatrix->SetElement(0, 3, 0.0);
+        rotationMatrix->SetElement(1, 3, 0.0);
+        rotationMatrix->SetElement(2, 3, 0.0);
+    }
 }
 
 // -------------------- initLookupTable --------------------

sdk/libraries/core/component/image/ImageComponent.h

@@ -216,6 +216,8 @@ public:
     virtual QWidget * getPropertyWidgetAt( unsigned int i, QWidget* parent = 0 );
     ///@}
 
+    const vtkSmartPointer<vtkMatrix4x4> getRotationMatrix() { return rotationMatrix; }
+
 protected:
 
     /**
@@ -223,7 +225,7 @@ protected:
      * @param anImageData The main vtkImageData of the volumic image.
      * @param copy Indicate if we do a vtk deep copy of these data or directly work on the one provided.
      * @param initialOrientation Initial image orientation
-     * @param initialTransformMatrix Initial image rotation (provided as a 3x3 matrix)
+     * @param initialTransformMatrix Initial image rotation (provided as a 4x4 matrix)
      */
     virtual void setImageData( vtkSmartPointer<vtkImageData> anImageData,
                                bool copy,
@@ -302,14 +304,20 @@ private:
     /// Initial image orientation
     ImageOrientationHelper::PossibleImageOrientations initialOrientation;
 
-    /// The transform to get the image back to the state it was at opening
-    /// This transform is the inverse of the transform of these 3 steps:
-    /// * the initial translation of the image (sometimes call Offset)
-    /// * the initial rotation of the image (sometimes call TransformMatrix)
-    /// * the initial image orientation
-    /// Mainly used to prepare the image for saving
-    vtkSmartPointer<vtkTransform> backToOriginalTransform;
+    /// The initial transform to the vtkImageData
+    /// @note This transform is equal to the initial image translation (Offset)
+    /// multiplies by the reorientation transform (transform -> RAI)
+    vtkSmartPointer<vtkTransform> initialImageDataTransform;
+
+    /// The initial frame of the image at opening
+    /// @note This transform is equal to the initial image translation (Offset)
+    /// multiplies by its rotation (TransformMatrix)
+    /// multiplies by the reorientation transform (transform -> RAI)
+    vtkSmartPointer<vtkTransform> initialFrameTransform;
 
+    /// The rotation matrix, that might have been altered by the user
+    /// Will be saved in header file information as TransformMatrix tag.
+    vtkSmartPointer<vtkMatrix4x4> rotationMatrix;
 };
 
 }