Skip to content
Snippets Groups Projects
Commit 3a8fabae authored by Mark Jenkinson's avatar Mark Jenkinson
Browse files

Moved orientation function

parent 66e63817
No related branches found
No related tags found
No related merge requests found
...@@ -1030,98 +1030,28 @@ float rms_deviation(const Matrix& affmat1, const Matrix& affmat2, ...@@ -1030,98 +1030,28 @@ float rms_deviation(const Matrix& affmat1, const Matrix& affmat2,
} }
Matrix mat44_to_newmat(mat44 inmat) // helper function - calls nifti, but with FSL default case
void get_axis_orientations(const Matrix& sform_mat, int sform_code,
const Matrix& qform_mat, int qform_code,
int& icode, int& jcode, int& kcode)
{ {
Matrix retmat(4,4); Matrix vox2mm(4,4);
for (int ii=0; ii<4; ii++) { if (sform_code!=NIFTI_XFORM_UNKNOWN) {
for (int jj=0; jj<4; jj++) { vox2mm = sform_mat;
retmat(ii+1,jj+1) = inmat.m[ii][jj]; } else if (qform_code!=NIFTI_XFORM_UNKNOWN) {
} vox2mm = qform_mat;
} } else {
return retmat; // ideally should be sampling_mat(), but for orientation it doesn't matter
} vox2mm = Identity(4);
vox2mm(1,1) = -vox2mm(1,1);
}
mat44 newmat_to_mat44(const Matrix& inmat) mat44 v2mm;
{ for (int ii=0; ii<4; ii++) { for (int jj=0; jj<4; jj++) {
mat44 retmat; v2mm.m[ii][jj] = vox2mm(ii+1,jj+1);
for (int ii=0; ii<4; ii++) { } }
for (int jj=0; jj<4; jj++) { mat44_to_orientation(v2mm,&icode,&jcode,&kcode);
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 // Added by MWW
......
...@@ -24,10 +24,11 @@ ...@@ -24,10 +24,11 @@
#include <string> #include <string>
#include <iterator> #include <iterator>
#include "fslio/fslio.h"
#include "config.h" #include "config.h"
#include "newmatap.h" #include "newmatap.h"
#include "kernel.h" #include "kernel.h"
#include "fslio/fslio.h"
//#pragma interface //#pragma interface
...@@ -157,27 +158,11 @@ namespace MISCMATHS { ...@@ -157,27 +158,11 @@ namespace MISCMATHS {
const ColumnVector& centre, const float rmax); const ColumnVector& centre, const float rmax);
float rms_deviation(const Matrix& affmat1, const Matrix& affmat2, float rms_deviation(const Matrix& affmat1, const Matrix& affmat2,
const float rmax=80.0); 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);
void get_axis_orientations(const Matrix& sform_mat, int sform_code,
const Matrix& qform_mat, int qform_code,
int& icode, int& jcode, int& kcode);
float median(const ColumnVector& x); float median(const ColumnVector& x);
void cart2sph(const ColumnVector& dir, float& th, float& ph);// cartesian to sperical polar coordinates void cart2sph(const ColumnVector& dir, float& th, float& ph);// cartesian to sperical polar coordinates
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment