src/gdcmOrientation.cxx

   1 /*=========================================================================
   2
   3   Program:   gdcm
   4   Module:    $RCSfile: gdcmOrientation.cxx,v $
   5   Language:  C++
   6   Date:      $Date: 2005/10/01 19:39:16 $
   7   Version:   $Revision: 1.15 $
   8
   9   Copyright (c) CREATIS (Centre de Recherche et d'Applications en Traitement de
  10   l'Image). All rights reserved. See Doc/License.txt or
  11   http://www.creatis.insa-lyon.fr/Public/Gdcm/License.html for details.
  12
  13      This software is distributed WITHOUT ANY WARRANTY; without even
  14      the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  15      PURPOSE.  See the above copyright notices for more information.
  16
  17 =========================================================================*/
  18
  19 #include "gdcmOrientation.h"
  20 #include "gdcmFile.h"
  21 #include "gdcmDebug.h"
  22 #include <math.h> // for sqrt
  23
  24 namespace gdcm
  25 {
  26 //--------------------------------------------------------------------
  27 //  THERALYS Algorithm to determine the most similar basic orientation
  28 //
  29 //  Transliterated from original Python code.
  30 //  Kept as close as possible to the original code
  31 //  in order to speed up any further modif of Python code :-(
  32 //-----------------------------------------------------------------------
  33
  34 /**
  35  * \brief  THERALYS' Algorithm to determine the most similar basic orientation
  36  *           (Axial, Coronal, Sagital) of the image
  37  * \note Should be run on the first gdcm::File of a 'coherent' Serie
  38  * @return orientation code
  39  *   #   0 : Not Applicable (neither 0020,0037 Image Orientation Patient
  40  *   #                       nor     0020,0032 Image Position           found)
  41  *   #   1 : Axial
  42  *   #  -1 : Axial invert
  43  *   #   2 : Coronal
  44  *   #  -2 : Coronal invert
  45  *   #   3 : Sagital
  46  *   #  -3 : Sagital invert
  47  *   #   4 : Heart Axial
  48  *   #  -4 : Heart Axial invert
  49  *   #   5 : Heart Coronal
  50  *   #  -5 : Heart Coronal invert
  51  *   #   6 : Heart Sagital
  52  *   #  -6 : Heart Sagital invert
  53  */
  54
  55 static const char  *OrientationTypeStrings[] = {
  56   "Not Applicable",
  57   "Axial",
  58   "Coronal",
  59   "Sagital",
  60   "Heart Axial",
  61   "Heart Coronal",
  62   "Heart Sagital",
  63   "Axial invert",
  64   "Coronal invert",
  65   "Sagital invert",
  66   "Heart Axial invert",
  67   "Heart Coronal invert",
  68   "Heart Sagital invert",
  69   NULL
  70 };
  71
  72 const char* Orientation::GetOrientationTypeString(OrientationType const o)
  73 {
  74   int k = (int)o;
  75   if (k < 0)
  76        k = -k + 6;
  77
  78   return OrientationTypeStrings[k];
  79 }
  80
  81 OrientationType Orientation::GetOrientationType( File *f )
  82 {
  83    float iop[6];
  84    bool succ = f->GetImageOrientationPatient( iop );
  85    if ( !succ )
  86    {
  87       gdcmErrorMacro( "No Image Orientation (0020,0037)/(0020,0032) found in the file, cannot proceed." )
  88       return NotApplicable;
  89    }
  90    vector3D ori1;
  91    vector3D ori2;
  92
  93    ori1.x = iop[0]; ori1.y = iop[1]; ori1.z = iop[2];
  94    ori2.x = iop[3]; ori2.y = iop[4]; ori2.z = iop[5];
  95
  96    // two perpendicular vectors describe one plane
  97    double dicPlane[6][2][3] =
  98    { {  { 1,   0,    0   },{ 0,      1,     0     }  }, // Axial
  99      {  { 1,   0,    0   },{ 0,      0,    -1     }  }, // Coronal
 100      {  { 0,   1,    0   },{ 0,      0,    -1     }  }, // Sagittal
 101      {  { 0.8, 0.5,  0.0 },{-0.1,    0.1 , -0.95  }  }, // Axial - HEART
 102      {  { 0.8, 0.5,  0.0 },{-0.6674, 0.687, 0.1794}  }, // Coronal - HEART
 103      {  {-0.1, 0.1, -0.95},{-0.6674, 0.687, 0.1794}  }  // Sagittal - HEART
 104    };
 105
 106    vector3D refA;
 107    vector3D refB;
 108    int i = 0;
 109    Res res;   // [ <result> , <memory of the last succes calcule> ]
 110    res.first = 0;
 111    res.second = 99999;
 112
 113    for (int numDicPlane=0; numDicPlane<6; numDicPlane++)
 114    {
 115        ++i;
 116        // refA=plane[0]
 117        refA.x = dicPlane[numDicPlane][0][0];
 118        refA.y = dicPlane[numDicPlane][0][1];
 119        refA.z = dicPlane[numDicPlane][0][2];
 120        // refB=plane[1]
 121        refB.x = dicPlane[numDicPlane][1][0];
 122        refB.y = dicPlane[numDicPlane][1][1];
 123        refB.z = dicPlane[numDicPlane][1][2];
 124        res=VerfCriterion(  i, CalculLikelyhood2Vec(refA,refB,ori1,ori2), res );
 125        res=VerfCriterion( -i, CalculLikelyhood2Vec(refB,refA,ori1,ori2), res );
 126    }
 127    gdcmAssertMacro( res.first <= 6 && res.first >= -6);
 128    return (OrientationType)res.first;
 129 }
 130
 131 Res
 132 Orientation::VerfCriterion(int typeCriterion, double criterionNew, Res const &in)
 133 {
 134    Res res;
 135    double type = in.first;
 136    double criterion = in.second;
 137    if (/*criterionNew < 0.1 && */criterionNew < criterion)
 138    {
 139       type      = typeCriterion;
 140       criterion = criterionNew;
 141    }
 142    res.first  = type;
 143    res.second = criterion;
 144    return res;
 145 }
 146
 147 inline double square_dist(vector3D const &v1, vector3D const &v2)
 148 {
 149   double res;
 150   res = (v1.x - v2.x)*(v1.x - v2.x) +
 151         (v1.y - v2.y)*(v1.y - v2.y) +
 152         (v1.z - v2.z)*(v1.z - v2.z);
 153   return res;
 154 }
 155
 156 //------------------------- Purpose : -----------------------------------
 157 //- This function determines the orientation similarity of two planes.
 158 //  Each plane is described by two vectors.
 159 //------------------------- Parameters : --------------------------------
 160 //- <refA>  : - type : vector 3D (double)
 161 //- <refB>  : - type : vector 3D (double)
 162 //            - Description of the first plane
 163 //- <ori1>  : - type : vector 3D (double)
 164 //- <ori2>  : - type : vector 3D (double)
 165 //            - Description of the second plane
 166 //------------------------- Return : ------------------------------------
 167 // double :   0 if the planes are perpendicular. While the difference of
 168 //            the orientation between the planes are big more enlarge is
 169 //            the criterion.
 170 //------------------------- Other : -------------------------------------
 171 // The calculus is based with vectors normalice
 172 double
 173 Orientation::CalculLikelyhood2Vec(vector3D const &refA, vector3D const &refB,
 174                                   vector3D const &ori1, vector3D const &ori2 )
 175 {
 176
 177    vector3D ori3 = ProductVectorial(ori1,ori2);
 178    vector3D refC = ProductVectorial(refA,refB);
 179    double res = square_dist(refC, ori3);
 180
 181    return sqrt(res);
 182 }
 183
 184 //------------------------- Purpose : -----------------------------------
 185 //- Calculus of the poduct vectorial between two vectors 3D
 186 //------------------------- Parameters : --------------------------------
 187 //- <vec1>  : - type : vector 3D (double)
 188 //- <vec2>  : - type : vector 3D (double)
 189 //------------------------- Return : ------------------------------------
 190 // (vec) :    - Vector 3D
 191 //------------------------- Other : -------------------------------------
 192 vector3D
 193 Orientation::ProductVectorial(vector3D const &vec1, vector3D const &vec2)
 194 {
 195    vector3D vec3;
 196    vec3.x =    vec1.y*vec2.z - vec1.z*vec2.y;
 197    vec3.y = -( vec1.x*vec2.z - vec1.z*vec2.x);
 198    vec3.z =    vec1.x*vec2.y - vec1.y*vec2.x;
 199
 200    return vec3;
 201 }
 202
 203
 204
 205 // ---------------------------------------------------------------------------
 206 // Here is the original Python code, kindly supplied by THERALYS
 207 //
 208 // C++ code doesn't give good results
 209 // --> FIXME
 210
 211 /*
 212
 213 def TypeOrientation(self,file0):
 214 """
 215 # ------------------------- Purpose : -----------------------------------
 216 # - This function compare the orientation of the given image and the
 217 #   basics orientations (Axial, Cornal, Sagital)
 218 # ------------------------- Parameters : --------------------------------
 219 # - <file0> : - type : string
 220 #             - The name of the first image file of the serie
 221 # ------------------------- Return : ------------------------------------
 222 #   1 :   Axial
 223 #  -1 :   Axial invert
 224 #   2 :   Coronal
 225 #  -2 :   Coronal invert
 226 #   3 :   Sagital
 227 #  -3 :   Sagital invert
 228 #   4 :   Heart Axial
 229 #  -4 :   Heart Axial invert
 230 #   5 :   Heart Coronal
 231 #  -5 :   Heart Coronal invert
 232 #   6 :   Heart Sagital
 233 #  -6 :   Heart Sagital invert
 234 #
 235    # ------------------------- Other : -------------------------------------
 236 # This method finds the most similar basic orientation.
 237 """
 238 try:
 239    toRead = gdcm.File(file0)
 240    ValDict = GetValuesDict(toRead)
 241    try:
 242       imageOrientation=ValDict["Image Orientation (Patient)"]
 243    except KeyError:
 244       imageOrientation=ValDict["Image Orientation"]
 245
 246    ori1=[float(split(imageOrientation,"\\")[0]),\
 247       float(split(imageOrientation,"\\")[1]),\
 248       float(split(imageOrientation,"\\")[2])]
 249    ori2=[float(split(imageOrientation,"\\")[3]),\
 250       float(split(imageOrientation,"\\")[4]),\
 251       float(split(imageOrientation,"\\")[5])]
 252
 253 ## two vectors perpendicular describe one plane
 254    dicPlane=[ [  [1,0,0],[0,1,0]   ],  ## Axial
 255             [  [1,0,0],[0,0,-1]  ],  ## Coronal
 256             [  [0,1,0],[0,0,-1]  ],  ## Sagittal
 257             [  [ 0.8 , 0.5 ,  0.0 ],[-0.1 , 0.1 , -0.95]        ],## Axial - HEART
 258             [  [ 0.8 , 0.5 ,  0.0 ],[-0.6674 , 0.687 , 0.1794]  ],## Coronal - HEART
 259             [  [-0.1 , 0.1 , -0.95],[-0.6674 , 0.687 , 0.1794]  ] ] ## Sagittal - HEART
 260
 261    i=0
 262    res=[0,99999]  ## [ <result> , <memory of the last succes calcule> ]
 263    for plane in dicPlane:
 264       i=i+1
 265       refA=plane[0]
 266       refB=plane[1]
 267       res=self.VerfCriterion(  i , self.CalculLikelyhood2Vec(refA,refB,ori1,ori2) , res )
 268       res=self.VerfCriterion( -i , self.CalculLikelyhood2Vec(refB,refA,ori1,ori2) , res )
 269    return res[0]
 270
 271    except KeyError:
 272    return 0
 273
 274
 275    def VerfCriterion(self,typeCriterion,criterionNew,res):
 276       type = res[0]
 277       criterion = res[1]
 278 #     if criterionNew<0.1 and criterionNew<criterion:
 279       if criterionNew<criterion:
 280          criterion=criterionNew
 281          type=typeCriterion
 282       return [ type , criterion ]
 283
 284
 285    def CalculLikelyhood2Vec(self,refA,refB,ori1,ori2):
 286 """
 287    # ------------------------- Purpose : -----------------------------------
 288    # - This function determine the orientation similarity of two planes.
 289    #   Each plane is described by two vector.
 290    # ------------------------- Parameters : --------------------------------
 291    # - <refA>  : - type : vector 3D (float)
 292    # - <refB>  : - type : vector 3D (float)
 293    #             - Description of the first plane
 294    # - <ori1>  : - type : vector 3D (float)
 295    # - <ori2>  : - type : vector 3D (float)
 296    #             - Description of the second plane
 297    # ------------------------- Return : ------------------------------------
 298    #  float :   0 if the planes are perpendicular.
 299    # While the difference of the orientation between the planes
 300    # are big more enlarge is
 301    # the criterion.
 302    # ------------------------- Other : -------------------------------------
 303    #  The calculus is based with vectors normalice
 304    """
 305
 306       ori3=self.ProductVectorial(ori1,ori2)
 307       refC=self.ProductVectorial(refA,refB)
 308       res=math.pow(refC[0]-ori3[0],2) + math.pow(refC[1]-ori3[1],2) + math.pow(refC[2]-ori3[2],2)
 309       return math.sqrt(res)
 310
 311    def ProductVectorial(self,vec1,vec2):
 312       """
 313       # ------------------------- Purpose : -----------------------------------
 314       # - Calculus of the poduct vectorial between two vectors 3D
 315       # ------------------------- Parameters : --------------------------------
 316       # - <vec1>  : - type : vector 3D (float)
 317       # - <vec2>  : - type : vector 3D (float)
 318       # ------------------------- Return : ------------------------------------
 319       #  (vec) :    - Vector 3D
 320       # ------------------------- Other : -------------------------------------
 321       """
 322       vec3=[0,0,0]
 323       vec3[0]=vec1[1]*vec2[2] - vec1[2]*vec2[1]
 324       vec3[1]=-( vec1[0]*vec2[2] - vec1[2]*vec2[0])
 325       vec3[2]=vec1[0]*vec2[1] - vec1[1]*vec2[0]
 326       return vec3
 327
 328    def GetValuesDict(image):
 329       """
 330       Returns a dictionnary containing values associated with Field Names
 331       dict["Dicom Field Name"]="Dicom field value"
 332       """
 333       val=image.GetFirstEntry()
 334       dic={}
 335       while(val):
 336          if isinstance(val,gdcm.ValEntryPtr):
 337             dic[val.GetName()]=val.GetValue()
 338          val=image.GetNextEntry()
 339       return dic
 340
 341 */
 342
 343
 344 // ------------------------------------------------------------------------
 345 /*
 346 2.2.2 Orientation of DICOM images
 347
 348
 349 http://www.dclunie.com/medical-image-faq/html/part2.html#DICOMOrientation
 350 says :
 351
 352 A question that is frequently asked in comp.protocols.dicom is how to determine
 353  which side of an image is which (e.g. left, right) and so on.
 354  The short answer is that for projection radiographs this is specified
 355  explicitly using the "Patient Orientation" attribute, and for cross-sectional
 356  images it needs to be derived from the "Image Orientation (Patient)" direction
 357  cosines. In the standard these are explained as follows:
 358
 359     * "C.7.6.1.1.1 Patient Orientation.
 360                 The Patient Orientation (0020,0020) relative to the image
 361                 plane shall be specified by two values that designate the
 362                 anatomical direction of the positive row axis (left to right)
 363                 and the positive column axis (top to bottom).
 364                 The first entry is the direction of the rows, given by the
 365                 direction of the last pixel in the first row from the first
 366                 pixel in that row.
 367                 The second entry is the direction of the columns, given by
 368                 the direction of the last pixel in the first column from the
 369                 first pixel in that column.
 370                 Anatomical direction shall be designated by the capital
 371                 letters: A (anterior), P (posterior), R (right),L (left),
 372                 H (head), F (foot).
 373                 Each value of the orientation attribute shall contain at
 374                 least one of these characters.
 375                 If refinements in the orientation descriptions are to be
 376                 specified, then they shall be designated by one or two
 377                 additional letters in each value.
 378                 Within each value, the letters shall be ordered with the
 379                 principal orientation designated in the first character."
 380
 381     * "C.7.6.2.1.1 Image Position And Image Orientation.
 382                 The "Image Position (Patient)" (0020,0032) specifies the x, y,
 383                 and z coordinates of the upper left hand corner of the image;
 384                 it is the center of the first voxel transmitted.
 385                 The "Image Orientation (Patient)" (0020,0037) specifies the
 386                 direction cosines of the first row and the first column
 387                 with respect to the patient.
 388                 These Attributes shall be provided as a pair.
 389                 Row value for the x, y, and z axes respectively followed by
 390                 the Column value for the x, y, and z axes respectively.
 391                 The direction of the axes is defined fully by the patient's
 392                 orientation.
 393                 The x-axis is increasing to the left hand side of the patient.
 394                 The y-axis is increasing to the posterior side of the patient.
 395                 The z-axis is increasing toward the head of the patient.
 396                 The patient based coordinate system is a right handed system,
 397                 i.e. the vector cross product of a unit vector along the
 398                 positive x-axis and a unit vector along the positive y-axis
 399                 is equal to a unit vector along the positive z-axis."
 400
 401 Some simple code to take one of the direction cosines (vectors) from the
 402 Image Orientation (Patient) attribute and generate strings equivalent to one
 403 of the values of Patient Orientation looks like this (noting that if the vector
 404 is not aligned exactly with one of the major axes, the resulting string will
 405 have multiple letters in as described under "refinements" in C.7.6.1.1.1):
 406
 407 */
 408
 409 /**
 410  * \brief computes the Patient Orientation relative to the image plane
 411  *          from the 'Image Orientation (Patient)'
 412  *          The first entry is the direction of the rows, given by the
 413  *          direction of the last pixel in the first row from the first
 414  *          pixel in that row.
 415  *          The second entry is the direction of the columns, given by
 416  *          the direction of the last pixel in the first column from the
 417  *          first pixel in that column.
 418  *          Anatomical direction is designated by the capital
 419  *          letters: A (anterior), P (posterior), R (right),L (left),
 420  *          H (head), F (foot).
 421  *          Refinements in the orientation descriptions are designated
 422  *          by one or two additional letters in each value.
 423  *          Use it when "Patient Orientation" (0020,0020) is not found
 424  * @return orientation string as "rawOrientation\columnsOrientation"
 425  */
 426 std::string Orientation::GetOrientation ( File *f )
 427 {
 428    float iop[6];
 429    if ( !f->GetImageOrientationPatient( iop ) )
 430    return GDCM_UNFOUND;
 431
 432    std::string orientation;
 433    orientation = GetSingleOrientation ( iop )
 434                + "\\"
 435                + GetSingleOrientation ( iop + 3 );
 436    return orientation;
 437 }
 438
 439
 440 std::string Orientation::GetSingleOrientation ( float *iop)
 441 {
 442    std::string orientation;
 443
 444    char orientationX = iop[0] < 0 ? 'R' : 'L';
 445    char orientationY = iop[1] < 0 ? 'A' : 'P';
 446    char orientationZ = iop[2] < 0 ? 'F' : 'H';
 447
 448    double absX = iop[0];
 449    if (absX < 0) absX = -absX;
 450       double absY = iop[1];
 451    if (absY < 0) absY = -absY;
 452       double absZ = iop[2];
 453    if (absZ < 0) absZ = -absZ;
 454
 455    for (int i=0; i<3; ++i)
 456    {
 457       if (absX>.0001 && absX>absY && absX>absZ)
 458       {
 459          orientation = orientation + orientationX;
 460          absX=0;
 461        }
 462        else if (absY>.0001 && absY>absX && absY>absZ)
 463        {
 464           orientation = orientation + orientationY;
 465           absY=0;
 466        }
 467        else if (absZ>.0001 && absZ>absX && absZ>absY)
 468        {
 469            orientation = orientation + orientationZ;
 470            absZ=0;
 471        }
 472        else
 473           break;
 474      }
 475    return orientation;
 476 }
 477
 478
 479 /*-------------------------------------------------------------------
 480
 481 Some more stuff, from XMedcon
 482
 483 ---> Personal remark from JPRx :
 484 --> patient_position (0x0018,0x5100) can be "HFS ", "FFS ", "HFP ", "FFP "
 485 --> the cosines may ahave any value -1< <+1, for MR images !
 486
 487 enum patient_slice_orientation_type
 488   {
 489     patient_slice_orientation_unknown = 0,
 490     supine_headfirst_transaxial,
 491     supine_headfirst_sagittal,
 492     supine_headfirst_coronal,
 493     supine_feetfirst_transaxial,
 494     supine_feetfirst_sagittal,
 495     supine_feetfirst_coronal,
 496     prone_headfirst_transaxial,
 497     prone_headfirst_sagittal,
 498     prone_headfirst_coronal,
 499     prone_feetfirst_transaxial,
 500     prone_feetfirst_sagittal,
 501     prone_feetfirst_coronal
 502   };
 503
 504 void GetImageOrientationPatient(gdcm::File &h,F32 image_orientation_patient[6])
 505 {
 506   h.GetImageOrientationPatient(image_orientation_patient);
 507 }
 508
 509 #if 0
 510 //
 511 // this is all completely cribbed from the xmedcon library, since
 512 // we're trying to do what it does, mostly.
 513 patient_slice_orientation_type
 514 GetPatSliceOrient(gdcm::File &h)
 515 {
 516   F32 image_orientation_patient[6];
 517
 518   // protected, do it the hard way
 519   //  h.GetImageOrientationPatient(image_orientation_patient);
 520   GetImageOrientationPatient(h,image_orientation_patient);
 521
 522   enum { headfirst, feetfirst } patient_orientation;
 523   enum { supine, prone } patient_rotation;
 524   enum { transaxial, sagittal, coronal } slice_orientation;
 525
 526   std::string patient_position = h.GetEntryByNumber(0x0018,0x5100);
 527   if(patient_position == "gdcm::Unfound")
 528     {
 529     patient_position = "HF";
 530     }
 531   if(patient_position.find("HF") != std::string::npos)
 532     {
 533     patient_orientation = headfirst;
 534     }
 535   else if(patient_position.find("FF") != std::string::npos)
 536     {
 537     patient_orientation = feetfirst;
 538     }
 539
 540   if(patient_position.find("S") != std::string::npos)
 541     {
 542     patient_rotation = supine;
 543     }
 544   else if(patient_position.find("P") != std::string::npos)
 545     {
 546     patient_rotation = prone;
 547     }
 548
 549   if((image_orientation_patient[0] == 1 || image_orientation_patient[0] == -1) &&
 550      (image_orientation_patient[4] == +1 || image_orientation_patient[4] == -1))
 551     {
 552     slice_orientation = transaxial;
 553     }
 554   else if((image_orientation_patient[1] == 1 || image_orientation_patient[1] == -1) &&
 555           (image_orientation_patient[5] == +1 || image_orientation_patient[5] == -1))
 556     {
 557     slice_orientation = sagittal;
 558     }
 559   else if((image_orientation_patient[0] == 1 || image_orientation_patient[0] == -1) &&
 560           (image_orientation_patient[5] == +1 || image_orientation_patient[5] == -1))
 561     {
 562     slice_orientation = coronal;
 563     }
 564   //
 565   // combine results
 566   patient_slice_orientation_type patient_slice_orientation =
 567     patient_slice_orientation_unknown;
 568   switch (patient_rotation)
 569     {
 570     case supine:
 571       switch (patient_orientation)
 572         {
 573         case headfirst:
 574           switch (slice_orientation)
 575             {
 576             case transaxial:
 577               patient_slice_orientation = supine_headfirst_transaxial;
 578               break;
 579             case sagittal:
 580               patient_slice_orientation = supine_headfirst_sagittal;
 581               break;
 582             case coronal:
 583               patient_slice_orientation = supine_headfirst_coronal;
 584               break;
 585             }
 586           break;
 587         case feetfirst:
 588           switch (slice_orientation)
 589             {
 590             case transaxial:
 591               patient_slice_orientation = supine_feetfirst_transaxial;
 592               break;
 593             case sagittal:
 594               patient_slice_orientation = supine_feetfirst_sagittal;
 595               break;
 596             case coronal:
 597               patient_slice_orientation = supine_feetfirst_coronal;
 598               break;
 599             }
 600           break;
 601         }
 602       break;
 603     case prone:
 604       switch (patient_orientation)
 605         {
 606         case headfirst:
 607           switch (slice_orientation)
 608             {
 609             case transaxial:
 610               patient_slice_orientation = prone_headfirst_transaxial;
 611               break;
 612             case sagittal:
 613               patient_slice_orientation = prone_headfirst_sagittal;
 614               break;
 615             case coronal:
 616               patient_slice_orientation = prone_headfirst_coronal;
 617               break;
 618             }
 619           break;
 620         case feetfirst:
 621           switch (slice_orientation)
 622             {
 623             case transaxial:
 624               patient_slice_orientation = prone_feetfirst_transaxial;
 625               break;
 626             case sagittal:
 627               patient_slice_orientation = prone_feetfirst_sagittal;
 628               break;
 629             case coronal:
 630               patient_slice_orientation = prone_feetfirst_coronal;
 631               break;
 632             }
 633           break;
 634         }
 635       break;
 636     }
 637   if(patient_slice_orientation != patient_slice_orientation_unknown)
 638     return patient_slice_orientation;
 639   //
 640   // this is what xmedcon does
 641   if ((image_orientation_patient[0] == +1)   &&
 642       (image_orientation_patient[4] == +1))
 643     patient_slice_orientation = supine_headfirst_transaxial;
 644   else if ((image_orientation_patient[0] == -1)   &&
 645            (image_orientation_patient[4] == +1))
 646     patient_slice_orientation = supine_feetfirst_transaxial;
 647   else if ((image_orientation_patient[0] == -1)   &&
 648            (image_orientation_patient[4] == -1))
 649     patient_slice_orientation = prone_headfirst_transaxial;
 650   else if ((image_orientation_patient[0] == +1)   &&
 651            (image_orientation_patient[4] == -1))
 652     patient_slice_orientation = prone_feetfirst_transaxial;
 653
 654   else if ((image_orientation_patient[1] == +1)   &&
 655            (image_orientation_patient[5] == -1))
 656     patient_slice_orientation = supine_headfirst_sagittal;
 657   else if ((image_orientation_patient[1] == +1)   &&
 658            (image_orientation_patient[5] == +1))
 659     patient_slice_orientation = supine_feetfirst_sagittal;
 660   else if ((image_orientation_patient[1] == -1)   &&
 661            (image_orientation_patient[5] == -1))
 662     patient_slice_orientation = prone_headfirst_sagittal;
 663   else if ((image_orientation_patient[1] == -1)   &&
 664            (image_orientation_patient[5] == +1))
 665     patient_slice_orientation = prone_feetfirst_sagittal;
 666
 667   else if ((image_orientation_patient[0] == +1)   &&
 668            (image_orientation_patient[5] == -1))
 669     patient_slice_orientation = supine_headfirst_coronal;
 670   else if ((image_orientation_patient[0] == -1)   &&
 671            (image_orientation_patient[5] == +1))
 672     patient_slice_orientation = supine_feetfirst_coronal;
 673   else if ((image_orientation_patient[0] == -1)   &&
 674            (image_orientation_patient[5] == -1))
 675     patient_slice_orientation = prone_headfirst_coronal;
 676   else if ((image_orientation_patient[0] == +1)   &&
 677            (image_orientation_patient[5] == +1))
 678     patient_slice_orientation = prone_feetfirst_coronal;
 679   return patient_slice_orientation;
 680 }
 681 #else
 682
 683 -------------------------------------------------------------------------*/
 684
 685 } // end namespace gdcm