Changes in coordinate lookup

miscmaths.cc

@@ -1030,6 +1030,99 @@ float rms_deviation(const Matrix& affmat1, const Matrix& affmat2,
 }
 
 
+Matrix mat44_to_newmat(mat44 inmat)
+{
+  Matrix retmat;
+  for (int ii=0; ii<4; ii++) {
+    for (int jj=0; jj<4; jj++) {
+      retmat(ii+1,jj+1) = inmat.m[ii][jj];
+    }
+  }
+  return retmat;
+}
+
+
+mat44 newmat_to_mat44(const Matrix& inmat)
+{
+  mat44 retmat;
+  for (int ii=0; ii<4; ii++) {
+    for (int jj=0; jj<4; jj++) {
+      retmat.m[ii][jj] = inmat(ii+1,jj+1);
+    }
+  }
+  return retmat;
+}
+
+
+int FslGetLeftRightOrder(int sform_code, const Matrix& sform_mat, 
+			 int qform_code, const Matrix& qform_mat)
+{
+  int retval;
+  // call the function within fslio
+  retval = FslGetLeftRightOrder2(sform_code,newmat_to_mat44(sform_mat),
+				 qform_code,newmat_to_mat44(qform_mat));
+  return retval;
+}
+
+short FslGetVox2mmMatrix(Matrix& vox2mm, 
+			 int sform_code, const Matrix& sform_mat, 
+			 int qform_code, const Matrix& qform_mat, 
+			 float dx, float dy, float dz)
+{
+  int retval;
+  mat44 vox2mm44;
+  // call the function within fslio
+  retval = FslGetVox2mmMatrix2(&vox2mm44,sform_code,newmat_to_mat44(sform_mat),
+			       qform_code,newmat_to_mat44(qform_mat),
+			       dx,dy,dz);
+  vox2mm = mat44_to_newmat(vox2mm44);
+  return retval;
+}
+
+
+
+Matrix Vox2FlirtCoord(int sform_code, const Matrix& sform_mat, 
+		      int qform_code, const Matrix& qform_mat, 
+		      float dx, float dy, float dz,
+		      int nx, int ny, int nz) 
+{
+  Matrix v2f(4,4);
+  Identity(v2f);
+  v2f(1,1)=dx;  v2f(2,2)=dy;  v2f(3,3)=dz;
+  if (FslGetLeftRightOrder(sform_code,sform_mat,qform_code,qform_mat)
+      == FSL_NEUROLOGICAL) {
+    Matrix swapx(4,4);
+    Identity(swapx);
+    swapx(1,1)=-1;
+    swapx(1,4)=nx-1;
+    v2f = v2f * swapx;
+  }
+  return v2f;
+}
+
+
+Matrix FslGetVox2VoxMatrix(const Matrix& flirt_in2ref,
+			   int sform_code_in, const Matrix& sform_mat_in, 
+			   int qform_code_in, const Matrix& qform_mat_in, 
+			   float dx_in, float dy_in, float dz_in,
+			   int nx_in, int ny_in, int nz_in,
+			   int sform_code_ref, const Matrix& sform_mat_ref, 
+			   int qform_code_ref, const Matrix& qform_mat_ref, 
+			   float dx_ref, float dy_ref, float dz_ref,
+			   int nx_ref, int ny_ref, int nz_ref)
+{
+  Matrix vox2flirt_in, vox2flirt_ref, vox2vox;
+  vox2flirt_in = Vox2FlirtCoord(sform_code_in,sform_mat_in,qform_code_in,
+				  qform_mat_in,dx_in,dy_in,dz_in,
+				  nx_in,ny_in,nz_in);
+  vox2flirt_ref = Vox2FlirtCoord(sform_code_ref,sform_mat_ref,qform_code_ref,
+				  qform_mat_ref,dx_ref,dy_ref,dz_ref,
+				  nx_ref,ny_ref,nz_ref);
+  vox2vox = vox2flirt_ref.i() * flirt_in2ref * vox2flirt_in;
+  return vox2vox;
+}
+
+
 // Added by MWW
 
 //  int getdiag(ColumnVector& diagvals, const Matrix& m)

miscmaths.h

@@ -27,6 +27,7 @@
 #include "config.h"
 #include "newmatap.h"
 #include "kernel.h"
+#include "fslio/fslio.h"
 
 //#pragma interface
 
@@ -156,6 +157,27 @@ namespace MISCMATHS {
 		      const ColumnVector& centre, const float rmax); 
   float rms_deviation(const Matrix& affmat1, const Matrix& affmat2, 
 		      const float rmax=80.0); 
+  // the following give left-right order and voxel2mm coordinate conversions
+  //  consistent with FSLView behaviour (NB: you cannot tell left-right from
+  //  the vox2mm matrix!)
+  int FslGetLeftRightOrder(int sform_code, const Matrix& sform_mat, 
+			   int qform_code, const Matrix& qform_mat);
+  short FslGetVox2mmMatrix(Matrix& vox2mm,
+			   int sform_code, const Matrix& sform_mat, 
+			   int qform_code, const Matrix& qform_mat, 
+			   float dx, float dy, float dz);
+  Matrix FslGetVox2VoxMatrix(const Matrix& flirt_in2ref,
+			     int sform_code_in, const Matrix& sform_mat_in, 
+			     int qform_code_in, const Matrix& qform_mat_in, 
+			     float dx_in, float dy_in, float dz_in,
+			     int nx_in, int ny_in, int nz_in,
+			     int sform_code_ref, const Matrix& sform_mat_ref, 
+			     int qform_code_ref, const Matrix& qform_mat_ref, 
+			     float dx_ref, float dy_ref, float dz_ref,
+			     int nx_ref, int ny_ref, int nz_ref);
+  mat44 newmat_to_mat44(const Matrix& inmat);
+  Matrix mat44_to_newmat(mat44 inmat);
+
   float median(const ColumnVector& x);
   
   void cart2sph(const ColumnVector& dir, float& th, float& ph);// cartesian to sperical polar coordinates