]> Creatis software - gdcm.git/blobdiff - src/gdcmFile.cxx
COMP: Try to get read of the float issues
[gdcm.git] / src / gdcmFile.cxx
index 214d8035162bb3cdb4828970230063a55e1cf02e..e21f6e0a81c76b5da44d9acf5c78db0776194fe4 100644 (file)
@@ -3,8 +3,8 @@
   Program:   gdcm
   Module:    $RCSfile: gdcmFile.cxx,v $
   Language:  C++
-  Date:      $Date: 2005/07/07 16:37:40 $
-  Version:   $Revision: 1.249 $
+  Date:      $Date: 2005/07/24 02:10:48 $
+  Version:   $Revision: 1.262 $
                                                                                 
   Copyright (c) CREATIS (Centre de Recherche et d'Applications en Traitement de
   l'Image). All rights reserved. See Doc/License.txt or
 //                                   as the Z coordinate, 
 // 0. for all the coordinates if nothing is found
 //
+// Image Position (Patient) (0020,0032) VM=3 What is it used for?
+// -->
+//  The attribute Patient Orientation (0020,0020) from the General Image Module 
+// is of type 2C and has the condition Required if image does not require 
+// Image Orientation (0020,0037) and Image Position (0020,0032). 
+// However, if the image does require the attributes 
+// - Image Orientation (Patient) (0020,0037), VM=6
+// - Image Position Patient (0020,0032), VM=3
+// then attribute Patient Orientation (0020,0020) should not be present
+//  in the images.
+//
+// Remember also :
+// Patient Position (0018,5100) values : HFP   = Head First-Prone
+//                                       HFS   = Head First-Supine
+//                                       HFDR  = Head First-Decubitus Right
+//                                       HFDL = Head First-Decubitus Left
+//                                       FFDR = Feet First-Decubitus Right
+//                                       FFDL = Feet First-Decubitus Left
+//                                       FFP  = Feet First-Prone
+//                                       FFS  = Feet First-Supine
+//                    can also find      SEMIERECT
+//                                       SUPINE
+// CS 2 Patient Orientation (0020 0020)
+//               When the coordinates of the image 
+//               are always present, this field is almost never used.
+//               Better we don't tust it too much ...
+//               Found Values are :      L\P
+//                                       L\FP
+//                                       P\F
+//                                       L\F
+//                                       P\FR
+//                                       R\F
+//
+// (0020|0037) [Image Orientation (Patient)] [1\0\0\0\1\0 ]
+
+                                      
 // ---------------------------------------------------------------
 //
 #include "gdcmFile.h"
 #include <vector>
 #include <stdio.h> //sscanf
 #include <stdlib.h> // for atoi
+#include <math.h> // for pow
 
 namespace gdcm 
 {
+
 //-----------------------------------------------------------------------------
 // Constructor / Destructor
 
@@ -60,18 +98,6 @@ File::File():
    NumPixel = 0x0010;
 }
 
-/**
- * \brief  Constructor 
- * @param  filename name of the file whose header we want to analyze
- */
-File::File( std::string const &filename )
-     :Document( )
-{    
-   RLEInfo  = new RLEFramesInfo;
-   JPEGInfo = new JPEGFragmentsInfo;
-
-   Load( filename ); // gdcm::Document is first Loaded, then the 'File part'
-}
 
 /**
  * \brief   Canonical destructor.
@@ -98,20 +124,6 @@ bool File::Load( )
 
     return DoTheLoadingJob( );   
 }
-/**
- * \brief   Loader. (DEPRECATED : not to break the API)
- * @param   fileName file to be open for parsing
- * @return false if file cannot be open or no swap info was found,
- *         or no tag was found.
- */
-bool File::Load( std::string const &fileName ) 
-{
-   SetFileName( fileName );
-   if ( ! this->Document::Load( ) )
-      return false;
-
-   return DoTheLoadingJob( );
-}
 
 /**
  * \brief   Does the Loading Job (internal use only)
@@ -724,7 +736,7 @@ float File::GetZOrigin()
   * @param iop adress of the (6)float array to receive values
   * @return cosines of image orientation patient
   */
-void File::GetImageOrientationPatient( float iop[6] )
+bool File::GetImageOrientationPatient( float iop[6] )
 {
    std::string strImOriPat;
    //iop is supposed to be float[6]
@@ -737,6 +749,7 @@ void File::GetImageOrientationPatient( float iop[6] )
           &iop[0], &iop[1], &iop[2], &iop[3], &iop[4], &iop[5]) != 6 )
       {
          gdcmWarningMacro( "Wrong Image Orientation Patient (0020,0037). Less than 6 values were found." );
+         return false;
       }
    }
    //For ACR-NEMA
@@ -747,8 +760,10 @@ void File::GetImageOrientationPatient( float iop[6] )
           &iop[0], &iop[1], &iop[2], &iop[3], &iop[4], &iop[5]) != 6 )
       {
          gdcmWarningMacro( "wrong Image Orientation Patient (0020,0035). Less than 6 values were found." );
+         return false;
       }
    }
+   return true;
 }
 
 /**
@@ -1300,9 +1315,11 @@ void File::AnonymizeNoLoad()
       if ( d == NULL)
          continue;
 
-      if ( dynamic_cast<BinEntry *>(d)
-        || dynamic_cast<SeqEntry *>(d) )
-         continue;
+         if ( dynamic_cast<SeqEntry *>(d) )
+         {
+            gdcmWarningMacro( "You cannot 'Anonymize a SeqEntry ");
+            continue;
+         }
 
       offset = d->GetOffset();
       lgth =   d->GetLength();
@@ -1323,6 +1340,7 @@ void File::AnonymizeNoLoad()
 /**
  * \brief anonymize a File (remove Patient's personal info passed with
  *        AddAnonymizeElement()
+ * \note You cannot Anonymize a BinEntry (to be fixed)
  */
 bool File::AnonymizeFile()
 {
@@ -1361,11 +1379,19 @@ bool File::AnonymizeFile()
          if ( d == NULL)
             continue;
 
-         if ( dynamic_cast<BinEntry *>(d)
-           || dynamic_cast<SeqEntry *>(d) )
+         if ( dynamic_cast<SeqEntry *>(d) )
+         {
+            gdcmWarningMacro( "You cannot 'Anonymize' a SeqEntry ");
             continue;
+         }
 
-         SetValEntry ((*it).Value, (*it).Group, (*it).Elem);
+         if ( dynamic_cast<BinEntry *>(d) )
+         {
+            gdcmWarningMacro( "To 'Anonymize' a BinEntry, better use AnonymizeNoLoad (FIXME) ");
+            continue;
+         }
+         else
+            SetValEntry ((*it).Value, (*it).Group, (*it).Elem);
       }
 }
 
@@ -1725,6 +1751,210 @@ void File::ReadAndSkipEncapsulatedBasicOffsetTable()
    }
 }
 
+// These are the deprecated method that one day should be removed (after the next release)
+
+#ifndef GDCM_LEGACY_REMOVE
+/**
+ * \brief  Constructor (DEPRECATED : temporaryly kept not to break the API)
+ * @param  filename name of the file whose header we want to analyze
+ * @deprecated do not use any longer
+ */
+File::File( std::string const &filename )
+     :Document( )
+{    
+   RLEInfo  = new RLEFramesInfo;
+   JPEGInfo = new JPEGFragmentsInfo;
+
+   SetFileName( filename );
+   Load( ); // gdcm::Document is first Loaded, then the 'File part'
+}
+
+/**
+ * \brief   Loader. (DEPRECATED :  temporaryly kept not to break the API)
+ * @param   fileName file to be open for parsing
+ * @return false if file cannot be open or no swap info was found,
+ *         or no tag was found.
+ * @deprecated Use the Load() [ + SetLoadMode() ] + SetFileName() functions instead
+ */
+bool File::Load( std::string const &fileName ) 
+{
+   GDCM_LEGACY_REPLACED_BODY(File::Load(std::string), "1.2",
+                             File::Load());
+   SetFileName( fileName );
+   if ( ! this->Document::Load( ) )
+      return false;
+
+   return DoTheLoadingJob( );
+}
+#endif
+
+// -----------------------------------------------------------------------------------------
+//  THERALYS Algorithm to determine the most similar basic orientation
+//
+//  Transliterated from original Python code.
+//  Kept as close as possible to the original code
+//  in order to speed up any further modif of Python code :-(
+// ------------------------------------------------------------------------------------------
+
+/**
+ * \brief  THERALYS' Algorithm to determine the most similar basic orientation
+ *           (Axial, Coronal, Sagital) of the image
+ * \note Should be run on the first gdcm::File of a 'coherent' Serie
+ * @return orientation code
+ * @return orientation code
+ *   #   0 :   Not Applicable (neither 0020,0037 Image Orientation Patient 
+ *   #                         nor     0020,0032Image Position     found )
+ *   #   1 :   Axial
+ *   #  -1 :   Axial invert
+ *   #   2 :   Coronal
+ *   #  -2 :   Coronal invert
+ *   #   3 :   Sagital
+ *   #  -3 :   Sagital invert
+ *   #   4 :   Heart Axial
+ *   #  -4 :   Heart Axial invert
+ *   #   5 :   Heart Coronal
+ *   #  -5 :   Heart Coronal invert
+ *   #   6 :   Heart Sagital
+ *   #  -6 :   Heart Sagital invert
+ */
+double File::TypeOrientation( )
+{
+   float iop[6];
+   bool succ = GetImageOrientationPatient( iop );
+   if ( !succ )
+   {
+      return 0.;
+   }
+
+   vector3D ori1;
+   vector3D ori2;
+
+   ori1.x = iop[0]; ori1.y = iop[1]; ori1.z = iop[2]; 
+   ori1.x = iop[3]; ori2.y = iop[4]; ori2.z = iop[5];
+
+   // two perpendicular vectors describe one plane
+   double dicPlane[6][2][3] =
+   { {  {1,    0,    0   },{0,       1,     0     }  },       // Axial
+     {  {1,    0,    0   },{0,       0,    -1     }  },       // Coronal
+     {  {0,    1,    0   },{0,       0,    -1     }  },       // Sagittal
+     {  { 0.8, 0.5,  0.0 },{-0.1,    0.1 , -0.95  }  },       // Axial - HEART
+     {  { 0.8, 0.5,  0.0 },{-0.6674, 0.687, 0.1794}  },       // Coronal - HEART
+     {  {-0.1, 0.1, -0.95},{-0.6674, 0.687, 0.1794}  }        // Sagittal - HEART
+   };
+
+   vector3D refA;
+   vector3D refB;
+   int i = 0;
+   Res res;   // [ <result> , <memory of the last succes calcule> ]
+   res.first = 0;
+   res.second = 99999;
+   for (int numDicPlane=0; numDicPlane<6; numDicPlane++)
+   {
+       ++i;
+       // refA=plane[0]
+       refA.x = dicPlane[numDicPlane][0][0]; 
+       refA.y = dicPlane[numDicPlane][0][1]; 
+       refA.z = dicPlane[numDicPlane][0][2];
+       // refB=plane[1]
+       refB.x = dicPlane[numDicPlane][1][0]; 
+       refB.y = dicPlane[numDicPlane][1][1]; 
+       refB.z = dicPlane[numDicPlane][1][2];
+       res=VerfCriterion(  i, CalculLikelyhood2Vec(refA,refB,ori1,ori2), res );
+       res=VerfCriterion( -i, CalculLikelyhood2Vec(refB,refA,ori1,ori2), res );
+   }
+   return res.first;
+/*
+//             i=0
+//             res=[0,99999]  ## [ <result> , <memory of the last succes calculus> ]
+//             for plane in dicPlane:
+//                 i=i+1
+//                 refA=plane[0]
+//                 refB=plane[1]
+//                 res=self.VerfCriterion(  i , self.CalculLikelyhood2Vec(refA,refB,ori1,ori2) , res )
+//                 res=self.VerfCriterion( -i , self.CalculLikelyhood2Vec(refB,refA,ori1,ori2) , res )
+//             return res[0]
+*/
+
+}
+
+Res 
+File::VerfCriterion(int typeCriterion, double criterionNew, Res const & in)
+{
+   Res res;
+   double criterion = in.second;
+   if (criterionNew < criterion)
+   {
+      res.first  = criterionNew;
+      res.second = typeCriterion;
+   }
+/*
+//   type = res[0]
+//   criterion = res[1]
+// #     if criterionNew<0.1 and criterionNew<criterion:
+//   if criterionNew<criterion:
+//      criterion=criterionNew
+//      type=typeCriterion
+//   return [ type , criterion ]
+*/
+   return res;
+} 
+
+inline double square_dist(vector3D const &v1, vector3D const & v2)
+{
+  double res;
+  res = (v1.x - v2.x)*(v1.x - v2.x) +
+        (v1.y - v2.y)*(v1.y - v2.y) +
+        (v1.z - v2.z)*(v1.z - v2.z);
+  return res;
+}
+
+//------------------------- Purpose : -----------------------------------
+//- This function determines the orientation similarity of two planes.
+//  Each plane is described by two vectors.
+//------------------------- Parameters : --------------------------------
+//- <refA>  : - type : vector 3D (double)
+//- <refB>  : - type : vector 3D (double)
+//            - Description of the first plane
+//- <ori1>  : - type : vector 3D (double)
+//- <ori2>  : - type : vector 3D (double)
+//            - Description of the second plane
+//------------------------- Return : ------------------------------------
+// double :   0 if the planes are perpendicular. While the difference of
+//            the orientation between the planes are big more enlarge is
+//            the criterion.
+//------------------------- Other : -------------------------------------
+// The calculus is based with vectors normalice
+double
+File::CalculLikelyhood2Vec(vector3D const & refA, vector3D const & refB, 
+                           vector3D const & ori1, vector3D const & ori2 )
+{
+
+   vector3D ori3 = ProductVectorial(ori1,ori2);
+   vector3D refC = ProductVectorial(refA,refB);
+   double res = square_dist(refC, ori3);
+
+   return sqrt(res);
+}
+
+//------------------------- Purpose : -----------------------------------
+//- Calculus of the poduct vectorial between two vectors 3D
+//------------------------- Parameters : --------------------------------
+//- <vec1>  : - type : vector 3D (double)
+//- <vec2>  : - type : vector 3D (double)
+//------------------------- Return : ------------------------------------
+// (vec) :    - Vector 3D
+//------------------------- Other : -------------------------------------
+vector3D
+File::ProductVectorial(vector3D const & vec1, vector3D const & vec2)
+{
+   vector3D vec3;
+   vec3.x =    vec1.y*vec2.z - vec1.z*vec2.y;
+   vec3.y = -( vec1.x*vec2.z - vec1.z*vec2.x);
+   vec3.z =    vec1.x*vec2.y - vec1.y*vec2.x;
+
+   return vec3;
+}
+
 //-----------------------------------------------------------------------------
 // Print