From f43ff3e533c512f8a35d10644ad7dc6d373f14b7 Mon Sep 17 00:00:00 2001
From: Tim Behrens <behrens@fmrib.ox.ac.uk>
Date: Mon, 23 Aug 2004 12:52:36 +0000
Subject: [PATCH] *** empty log message ***
---
Makefile | 7 ++++-
make_dyadic_vectors.cc | 66 ++++++++++++++++++++++++++++++++++++++++++
2 files changed, 72 insertions(+), 1 deletion(-)
create mode 100644 make_dyadic_vectors.cc
diff --git a/Makefile b/Makefile
index b5da435..e4c633c 100644
--- a/Makefile
+++ b/Makefile
@@ -18,6 +18,7 @@ ROM=reord_OM
SAUS=sausages
DIFF_PVM=diff_pvm
RV=replacevols
+MDV=make_dyadic_vectors
DTIFITOBJS=dtifit.o dtifitOptions.o
PTOBJS=probtrack.o probtrackOptions.o pt_alltracts.o pt_matrix.o pt_seeds_to_targets.o pt_simple.o pt_twomasks.o pt_matrix_mesh.o
@@ -28,9 +29,10 @@ ROMOBJS=reord_OM.o
SAUSOBJS=sausages.o
DIFF_PVMOBJS=diff_pvm.o diff_pvmoptions.o
RVOBJS=replacevols.o
+MDVOBJS=make_dyadic_vectors.o
SCRIPTS = eddy_correct
-XFILES = dtifit probtrack find_the_biggest medianfilter diff_pvm
+XFILES = dtifit probtrack find_the_biggest medianfilter diff_pvm make_dyadic_vectors
RUNTCLS = Fdt
all: ${XFILES} reord_OM sausages
@@ -62,6 +64,9 @@ ${DIFF_PVM}: ${DIFF_PVMOBJS}
${RV}: ${RVOBJS}
${CXX} ${CXXFLAGS} ${LDFLAGS} -o $@ ${RVOBJS} ${DLIBS}
+${MDV}: ${MDVOBJS}
+ ${CXX} ${CXXFLAGS} ${LDFLAGS} -o $@ ${MDVOBJS} ${DLIBS}
+
diff --git a/make_dyadic_vectors.cc b/make_dyadic_vectors.cc
new file mode 100644
index 0000000..ad62ede
--- /dev/null
+++ b/make_dyadic_vectors.cc
@@ -0,0 +1,66 @@
+#include <iostream>
+#include <fstream>
+#include "newimage/newimageall.h"
+#include <vector>
+using namespace std;
+using namespace NEWIMAGE;
+
+int main ( int argc, char **argv ){
+ if(argc<4){
+ cerr<<"usage: make_dyadic_vectors <theta_vol> <phi_vol> <output>"<<endl;
+ exit(0);
+ }
+ volume4D<float> ths,phs;
+ volumeinfo tempinfo;
+ read_volume4D(ths,argv[1],tempinfo);
+ read_volume4D(phs,argv[2]);
+ volume4D<float> dyadic_vecs(ths.xsize(),ths.ysize(),ths.zsize(),3);
+ dyadic_vecs=0;
+ SymmetricMatrix dyad(3);dyad=0;
+ ColumnVector dir(3);
+
+ DiagonalMatrix dyad_D; //eigenvalues
+ Matrix dyad_V; //eigenvectors
+
+ for(int k=ths.minz();k<=ths.maxz();k++){
+ for(int j=ths.miny();j<=ths.maxy();j++){
+ for(int i=ths.minx();i<=ths.maxx();i++){
+ dyad=0;
+ for(int s=ths.mint();s<=ths.maxt();s++){
+ float th=ths(i,j,k,s);
+ float ph=phs(i,j,k,s);
+ dir(1)=sin(th)*cos(ph);
+ dir(2)=sin(th)*sin(ph);
+ dir(3)=cos(th);
+ dyad << dyad+dir*dir.t();
+ }
+
+ EigenValues(dyad,dyad_D,dyad_V);
+ int maxeig;
+ if(dyad_D(1)>dyad_D(2)){
+ if(dyad_D(1)>dyad_D(3)) maxeig=1;
+ else maxeig=3;
+ }
+ else{
+ if(dyad_D(2)>dyad_D(3)) maxeig=2;
+ else maxeig=3;
+ }
+ dyadic_vecs(i,j,k,0)=dyad_V(1,maxeig);
+ dyadic_vecs(i,j,k,1)=dyad_V(2,maxeig);
+ dyadic_vecs(i,j,k,2)=dyad_V(3,maxeig);
+ }
+
+ }
+ }
+
+ save_volume4D(dyadic_vecs,argv[3],tempinfo);
+}
+
+
+
+
+
+
+
+
+
--
GitLab