Program: gdcm
Module: $RCSfile: gdcmOrientation.cxx,v $
Language: C++
- Date: $Date: 2005/07/25 03:44:35 $
- Version: $Revision: 1.3 $
+ Date: $Date: 2005/09/19 09:48:27 $
+ Version: $Revision: 1.8 $
Copyright (c) CREATIS (Centre de Recherche et d'Applications en Traitement de
l'Image). All rights reserved. See Doc/License.txt or
* (Axial, Coronal, Sagital) of the image
* \note Should be run on the first gdcm::File of a 'coherent' Serie
* @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
+ * # 0 : Not Applicable (neither 0020,0037 Image Orientation Patient
+ * # nor 0020,0032 Image 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 Orientation::TypeOrientation( File *f )
{
bool succ = f->GetImageOrientationPatient( iop );
if ( !succ )
{
- gdcmErrorMacro( "No Image Orientation (0020,0037) was found in the file, cannot proceed." )
+ gdcmErrorMacro( "No Image Orientation (0020,0037) found in the file, cannot proceed." )
return 0;
}
-
+/*
+std::cout << " iop : ";
+for(int i=0;i<6;i++)
+ std::cout << iop[i] << " ";
+std::cout << std::endl;
+*/
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];
+ ori2.x = iop[3]; ori2.y = iop[4]; ori2.z = iop[5];
// two perpendicular vectors describe one plane
double dicPlane[6][2][3] =
}
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]
+// 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
-Orientation::VerfCriterion(int typeCriterion, double criterionNew, Res const & in)
+Orientation::VerfCriterion(int typeCriterion, double criterionNew, Res const &in)
{
Res res;
double criterion = in.second;
if (criterionNew < criterion)
{
- res.first = criterionNew;
- res.second = typeCriterion;
+ res.first = typeCriterion;;
+ res.second = criterionNew;
}
/*
// type = res[0]
return res;
}
-inline double square_dist(vector3D const &v1, vector3D const & v2)
+inline double square_dist(vector3D const &v1, vector3D const &v2)
{
double res;
res = (v1.x - v2.x)*(v1.x - v2.x) +
//------------------------- Other : -------------------------------------
// The calculus is based with vectors normalice
double
-Orientation::CalculLikelyhood2Vec(vector3D const & refA, vector3D const & refB,
- vector3D const & ori1, vector3D const & ori2 )
+Orientation::CalculLikelyhood2Vec(vector3D const &refA, vector3D const &refB,
+ vector3D const &ori1, vector3D const &ori2 )
{
vector3D ori3 = ProductVectorial(ori1,ori2);
} // end namespace gdcm
+
+
+
+// ---------------------------------------------------------------------------
+// Here is the original Python code, kindly supplied by THERALYS
+//
+// C++ code doesn't give good results
+// --> FIXME
+
+/*
+
+def TypeOrientation(self,file0):
+"""
+# ------------------------- Purpose : -----------------------------------
+# - This function compare the orientation of the given image and the
+# basics orientations (Axial, Cornal, Sagital)
+# ------------------------- Parameters : --------------------------------
+# - <file0> : - type : string
+# - The name of the first image file of the serie
+# ------------------------- Return : ------------------------------------
+# 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
+#
+ # ------------------------- Other : -------------------------------------
+# This method finds the most similar basic orientation.
+"""
+try:
+ toRead = gdcm.File(file0)
+ ValDict = GetValuesDict(toRead)
+ try:
+ imageOrientation=ValDict["Image Orientation (Patient)"]
+ except KeyError:
+ imageOrientation=ValDict["Image Orientation"]
+
+ ori1=[float(split(imageOrientation,"\\")[0]),\
+ float(split(imageOrientation,"\\")[1]),\
+ float(split(imageOrientation,"\\")[2])]
+ ori2=[float(split(imageOrientation,"\\")[3]),\
+ float(split(imageOrientation,"\\")[4]),\
+ float(split(imageOrientation,"\\")[5])]
+
+## two vectors perpendicular describe one plane
+ dicPlane=[ [ [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
+
+ i=0
+ res=[0,99999] ## [ <result> , <memory of the last succes calcule> ]
+ 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]
+
+ except KeyError:
+ return 0
+
+
+ def VerfCriterion(self,typeCriterion,criterionNew,res):
+ type = res[0]
+ criterion = res[1]
+# if criterionNew<0.1 and criterionNew<criterion:
+ if criterionNew<criterion:
+ criterion=criterionNew
+ type=typeCriterion
+ return [ type , criterion ]
+
+
+ def CalculLikelyhood2Vec(self,refA,refB,ori1,ori2):
+"""
+ # ------------------------- Purpose : -----------------------------------
+ # - This function determine the orientation similarity of two planes.
+ # Each plane is described by two vector.
+ # ------------------------- Parameters : --------------------------------
+ # - <refA> : - type : vector 3D (float)
+ # - <refB> : - type : vector 3D (float)
+ # - Description of the first plane
+ # - <ori1> : - type : vector 3D (float)
+ # - <ori2> : - type : vector 3D (float)
+ # - Description of the second plane
+ # ------------------------- Return : ------------------------------------
+ # float : 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
+ """
+
+ ori3=self.ProductVectorial(ori1,ori2)
+ refC=self.ProductVectorial(refA,refB)
+ res=math.pow(refC[0]-ori3[0],2) + math.pow(refC[1]-ori3[1],2) + math.pow(refC[2]-ori3[2],2)
+ return math.sqrt(res)
+
+ def ProductVectorial(self,vec1,vec2):
+ """
+ # ------------------------- Purpose : -----------------------------------
+ # - Calculus of the poduct vectorial between two vectors 3D
+ # ------------------------- Parameters : --------------------------------
+ # - <vec1> : - type : vector 3D (float)
+ # - <vec2> : - type : vector 3D (float)
+ # ------------------------- Return : ------------------------------------
+ # (vec) : - Vector 3D
+ # ------------------------- Other : -------------------------------------
+ """
+ vec3=[0,0,0]
+ vec3[0]=vec1[1]*vec2[2] - vec1[2]*vec2[1]
+ vec3[1]=-( vec1[0]*vec2[2] - vec1[2]*vec2[0])
+ vec3[2]=vec1[0]*vec2[1] - vec1[1]*vec2[0]
+ return vec3
+
+ def GetValuesDict(image):
+ """
+ Returns a dictionnary containing values associated with Field Names
+ dict["Dicom Field Name"]="Dicom field value"
+ """
+ val=image.GetFirstEntry()
+ dic={}
+ while(val):
+ if isinstance(val,gdcm.ValEntryPtr):
+ dic[val.GetName()]=val.GetValue()
+ val=image.GetNextEntry()
+ return dic
+
+*/