From 0baf4c8ba25e25dc9aa9298406f575f92777ad5e Mon Sep 17 00:00:00 2001
From: Tim Behrens <behrens@fmrib.ox.ac.uk>
Date: Fri, 15 Jul 2005 13:21:07 +0000
Subject: [PATCH] *** empty log message ***
---
probtrackx.cc | 1 -
probtrackx.h | 1 -
ptx_seedmask.cc | 1 -
ptx_seedmask.h | 2 +-
ptx_simple.cc | 1 -
ptx_simple.h | 2 +-
ptx_twomasks.cc | 1 -
ptx_twomasks.h | 2 +-
streamlines.cc | 2 +-
streamlines.h | 8 +-
tractvols.h | 197 +++++++++++++++++++------------------------
tractvolsx.h | 217 ++++++++++++++++++++++++++++++++++++++++++++++++
12 files changed, 309 insertions(+), 126 deletions(-)
create mode 100644 tractvolsx.h
diff --git a/probtrackx.cc b/probtrackx.cc
index 4c9cf19..03e7bbe 100644
--- a/probtrackx.cc
+++ b/probtrackx.cc
@@ -15,7 +15,6 @@ using namespace NEWIMAGE;
using namespace TRACT;
using namespace Utilities;
using namespace PARTICLE;
-using namespace TRACTVOLS;
using namespace mesh;
//using namespace NEWMAT;
//////////////////////////
diff --git a/probtrackx.h b/probtrackx.h
index 1e2fb56..2f52f3e 100644
--- a/probtrackx.h
+++ b/probtrackx.h
@@ -4,7 +4,6 @@
#include "probtrackxOptions.h"
#include "particle.h"
-#include "tractvols.h"
#include "ptx_simple.h"
#include "ptx_seedmask.h"
#include "ptx_twomasks.h"
diff --git a/ptx_seedmask.cc b/ptx_seedmask.cc
index 9539738..b771812 100644
--- a/ptx_seedmask.cc
+++ b/ptx_seedmask.cc
@@ -10,7 +10,6 @@ using namespace NEWIMAGE;
using namespace TRACT;
using namespace Utilities;
using namespace PARTICLE;
-using namespace TRACTVOLS;
using namespace mesh;
diff --git a/ptx_seedmask.h b/ptx_seedmask.h
index dbc4968..2c284c4 100644
--- a/ptx_seedmask.h
+++ b/ptx_seedmask.h
@@ -8,7 +8,7 @@
#include "meshclass/meshclass.h"
#include "probtrackxOptions.h"
#include "particle.h"
-#include "tractvols.h"
+
diff --git a/ptx_simple.cc b/ptx_simple.cc
index fd47364..5926273 100644
--- a/ptx_simple.cc
+++ b/ptx_simple.cc
@@ -10,7 +10,6 @@ using namespace NEWIMAGE;
using namespace TRACT;
using namespace Utilities;
using namespace PARTICLE;
-using namespace TRACTVOLS;
using namespace mesh;
diff --git a/ptx_simple.h b/ptx_simple.h
index eb57c1f..c9417e7 100644
--- a/ptx_simple.h
+++ b/ptx_simple.h
@@ -8,7 +8,7 @@
#include "meshclass/meshclass.h"
#include "probtrackxOptions.h"
#include "particle.h"
-#include "tractvols.h"
+
diff --git a/ptx_twomasks.cc b/ptx_twomasks.cc
index 6010297..68818b0 100644
--- a/ptx_twomasks.cc
+++ b/ptx_twomasks.cc
@@ -10,7 +10,6 @@ using namespace NEWIMAGE;
using namespace TRACT;
using namespace Utilities;
using namespace PARTICLE;
-using namespace TRACTVOLS;
using namespace mesh;
diff --git a/ptx_twomasks.h b/ptx_twomasks.h
index 40fce4a..1c767a2 100644
--- a/ptx_twomasks.h
+++ b/ptx_twomasks.h
@@ -8,6 +8,6 @@
#include "meshclass/meshclass.h"
#include "probtrackxOptions.h"
#include "particle.h"
-#include "tractvols.h"
+
void twomasks();
diff --git a/streamlines.cc b/streamlines.cc
index aa03ccc..a2691b2 100644
--- a/streamlines.cc
+++ b/streamlines.cc
@@ -73,7 +73,7 @@ namespace TRACT{
bool Streamliner::streamline(const float& x_init,const float& y_init,const float& z_init, const ColumnVector& dim_seeds,const int& fibst){
- //fibst tells tractvols which fibre to start with if there are more than one..
+ //fibst tells tractvolsx which fibre to start with if there are more than one..
//x_init etc. are in seed space...
vols.reset(fibst);
m_x_s_init=x_init; //seed x position in voxels
diff --git a/streamlines.h b/streamlines.h
index c7bf349..9fb8a75 100644
--- a/streamlines.h
+++ b/streamlines.h
@@ -4,12 +4,12 @@
#include "meshclass/meshclass.h"
#include "probtrackxOptions.h"
#include "particle.h"
-#include "tractvols.h"
+#include "tractvolsx.h"
using namespace std;
using namespace NEWIMAGE;
using namespace Utilities;
-using namespace TRACTVOLS;
+using namespace TRACTVOLSX;
using namespace mesh;
using namespace PARTICLE;
@@ -18,7 +18,7 @@ namespace TRACT{
class Streamliner{
- //Everything in DTI space is done INSIDE this class and lower level classes (particle and tractvols)
+ //Everything in DTI space is done INSIDE this class and lower level classes (particle and tractvolsx)
//This class communicates with higher level classes in Seed voxels.
//
probtrackxOptions& opts;
@@ -33,7 +33,7 @@ namespace TRACT{
vector<bool> m_passed_flags;
vector<bool> m_own_waymasks;
Matrix m_Seeds_to_DTI;
- TractVols vols;
+ Tractvolsx vols;
float m_lcrat;
float m_x_s_init;
float m_y_s_init;
diff --git a/tractvols.h b/tractvols.h
index b7a9b27..d743274 100644
--- a/tractvols.h
+++ b/tractvols.h
@@ -1,10 +1,70 @@
-/* tractvolsX.h
+/* tractvols.h
Tim Behrens, FMRIB Image Analysis Group
Copyright (C) 2004 University of Oxford */
-/* CCOPYRIGHT */
+/* Part of FSL - FMRIB's Software Library
+ http://www.fmrib.ox.ac.uk/fsl
+ fsl@fmrib.ox.ac.uk
+
+ Developed at FMRIB (Oxford Centre for Functional Magnetic Resonance
+ Imaging of the Brain), Department of Clinical Neurology, Oxford
+ University, Oxford, UK
+
+
+ LICENCE
+
+ FMRIB Software Library, Release 3.2beta (c) 2004, The University of
+ Oxford (the "Software")
+
+ The Software remains the property of the University of Oxford ("the
+ University").
+
+ The Software is distributed "AS IS" under this Licence solely for
+ non-commercial use in the hope that it will be useful, but in order
+ that the University as a charitable foundation protects its assets for
+ the benefit of its educational and research purposes, the University
+ makes clear that no condition is made or to be implied, nor is any
+ warranty given or to be implied, as to the accuracy of the Software,
+ or that it will be suitable for any particular purpose or for use
+ under any specific conditions. Furthermore, the University disclaims
+ all responsibility for the use which is made of the Software. It
+ further disclaims any liability for the outcomes arising from using
+ the Software.
+
+ The Licensee agrees to indemnify the University and hold the
+ University harmless from and against any and all claims, damages and
+ liabilities asserted by third parties (including claims for
+ negligence) which arise directly or indirectly from the use of the
+ Software or the sale of any products based on the Software.
+
+ No part of the Software may be reproduced, modified, transmitted or
+ transferred in any form or by any means, electronic or mechanical,
+ without the express permission of the University. The permission of
+ the University is not required if the said reproduction, modification,
+ transmission or transference is done without financial return, the
+ conditions of this Licence are imposed upon the receiver of the
+ product, and all original and amended source code is included in any
+ transmitted product. You may be held legally responsible for any
+ copyright infringement that is caused or encouraged by your failure to
+ abide by these terms and conditions.
+
+ You are not permitted under this Licence to use this Software
+ commercially. Use for which any financial return is received shall be
+ defined as commercial use, and includes (1) integration of all or part
+ of the source code or the Software into a product for sale or license
+ by or on behalf of Licensee to third parties or (2) use of the
+ Software or any derivative of it for research with the final aim of
+ developing software products for sale or license to a third party or
+ (3) use of the Software or any derivative of it for research with the
+ final aim of developing non-software products for sale or license to a
+ third party, or (4) use of the Software to provide any service to an
+ external organisation for which payment is received. If you are
+ interested in using the Software commercially, please contact Isis
+ Innovation Limited ("Isis"), the technology transfer company of the
+ University, to negotiate a licence. Contact details are:
+ innovation@isis.ox.ac.uk quoting reference DE/1112. */
#ifndef __TRACTVOLS_H_
#define __TRACTVOLS_H_
@@ -16,93 +76,36 @@
#include "newimage/newimageall.h"
#include <iostream>
#include "stdlib.h"
-#include "probtrackxOptions.h"
+
using namespace std;
using namespace NEWIMAGE;
-using namespace TRACT;
namespace TRACTVOLS{
class TractVols
{
private:
- probtrackxOptions& opts;
- vector<volume4D<float>* > thsamples;
- vector<volume4D<float>* > phsamples;
- vector<volume4D<float>* > fsamples;
- bool init_sample;
- int fibst;
+ volume4D<float> thsamples;
+ volume4D<float> phsamples;
+ volume4D<float> fsamples;
bool usef;
public:
//constructors::
- TractVols(const bool& usefin=false):opts(probtrackxOptions::getInstance()),init_sample(true),fibst(1),usef(usefin){}
- TractVols():opts(probtrackxOptions::getInstance()){}
- ~TractVols(){
- for(unsigned int m=0;m<thsamples.size();m++)
- delete thsamples[m]; //ask flitney, do you just delete the ptr??
- for(unsigned int m=0;m<phsamples.size();m++)
- delete phsamples[m];
- for(unsigned int m=0;m<fsamples.size();m++)
- delete fsamples[m];
- }
+ TractVols(const bool& usefin=false):usef(usefin){}
+ ~TractVols(){}
-
- void reset(const int& fibst_in){
- init_sample=true;
- fibst=fibst_in;
- }
//Initialise
void initialise(const string& basename){
-
-
- if(fsl_imageexists(basename+"_thsamples")){
- volume4D<float> *tmpthptr= new volume4D<float>;
- volume4D<float> *tmpphptr= new volume4D<float>;
- volume4D<float> *tmpfptr= new volume4D<float>;
- cout<<"1"<<endl;
- read_volume4D(*tmpthptr,basename+"_thsamples");
- cout<<"2"<<endl;
- thsamples.push_back(tmpthptr);
- cout<<"3"<<endl;
- read_volume4D(*tmpphptr,basename+"_phsamples");
- cout<<"4"<<endl;
- phsamples.push_back(tmpphptr);
- cout<<"5"<<endl;
- if(usef){
- read_volume4D(*tmpfptr,basename+"_fsamples");
- fsamples.push_back(tmpfptr);
- }
- cout<<"6"<<endl;
- }
- else{
- int fib=1;
- bool fib_existed=true;
- while(fib_existed){
- if(fsl_imageexists(basename+"_th"+num2str(fib)+"samples")){
- volume4D<float> *tmpthptr= new volume4D<float>;
- volume4D<float> *tmpphptr= new volume4D<float>;
- volume4D<float> *tmpfptr= new volume4D<float>;
- cout<<fib<<"_1"<<endl;
- read_volume4D(*tmpthptr,basename+"_th"+num2str(fib)+"samples");
- thsamples.push_back(tmpthptr);
- cout<<fib<<"_2"<<endl;
- read_volume4D(*tmpphptr,basename+"_ph"+num2str(fib)+"samples");
- phsamples.push_back(tmpphptr);
- cout<<fib<<"_3"<<endl;
- read_volume4D(*tmpfptr,basename+"_f"+num2str(fib)+"samples");
- fsamples.push_back(tmpfptr);
- fib++;
- }
- else{
- fib_existed=false;
- }
- }
-
- }
- cout<<"7"<<endl;
+ read_volume4D(thsamples,basename+"_thsamples");
+ read_volume4D(phsamples,basename+"_phsamples");
+ if(usef)
+ read_volume4D(fsamples,basename+"_fsamples");
}
-
- ColumnVector sample(const float& x,const float& y,const float&z,const float& r_x,const float& r_y,const float& r_z){
+
+ ColumnVector sample(const float& x,const float& y,const float&z){
+ // int r_x=(int) round(x);
+ // int r_y=(int) round(y);
+ // int r_z=(int) round(z);
////////Probabilistic interpolation
int cx =(int) ceil(x),fx=(int) floor(x);
@@ -150,49 +153,17 @@ namespace TRACTVOLS{
float samp=rand(); samp/=RAND_MAX;
- samp=round(samp*((*thsamples[0]).tsize()-1));
- float theta=0,phi=0;
- float dotmax=0,dottmp=0;
- int fibind=0;
- if(thsamples.size()>1){//more than 1 fibre
- if(init_sample){//go for the specified option on the first jump
- theta=(*thsamples[fibst])(int(newx),int(newy),int(newz),int(samp));
- phi=(*phsamples[fibst])(int(newx),int(newy),int(newz),int(samp));
- init_sample=false;
- }
- else{
- for(unsigned int fib=0;fib<thsamples.size();fib++){
- if((*fsamples[fib])(int(newx),int(newy),int(newz),int(samp))>opts.fibthresh.value()){
- float phtmp=(*phsamples[fib])(int(newx),int(newy),int(newz),int(samp));
- float thtmp=(*thsamples[fib])(int(newx),int(newy),int(newz),int(samp));
- dottmp=fabs(sin(thtmp)*cos(phtmp)*r_x + sin(thtmp)*sin(phtmp)*r_y + cos(thtmp)*r_z);
- if(dottmp>dotmax){
- dotmax=dottmp;
- theta=thtmp;
- phi=phtmp;
- fibind=fib;
- }
+ samp=round(samp*(thsamples.tsize()-1));
+ //float phi = phsamples(r_x,r_y,r_z,samp);
+ //float theta = thsamples(r_x,r_y,r_z,samp);
- }
-
- }
-
- if(dotmax==0){
- theta=(*thsamples[0])(int(newx),int(newy),int(newz),int(samp));
- phi=(*phsamples[0])(int(newx),int(newy),int(newz),int(samp));
- }
- }
- }
- else{
- theta=(*thsamples[0])(int(newx),int(newy),int(newz),int(samp));
- phi=(*phsamples[0])(int(newx),int(newy),int(newz),int(samp));
- }
-
+ float phi = phsamples(int(newx),int(newy),int(newz),int(samp));
+ float theta = thsamples(int(newx),int(newy),int(newz),int(samp));
float f;
if(usef){
- f = (*fsamples[fibind])(int(newx),int(newy),int(newz),int(samp));
+ f = fsamples(int(newx),int(newy),int(newz),int(samp));
}
else
f=1;
@@ -205,7 +176,7 @@ namespace TRACTVOLS{
ColumnVector dimensions(){
ColumnVector dims(3);
- dims << (*thsamples[0]).xdim() <<(*thsamples[0]).ydim() << (*thsamples[0]).zdim();
+ dims << thsamples.xdim() <<thsamples.ydim() << thsamples.zdim();
return dims;
}
};
diff --git a/tractvolsx.h b/tractvolsx.h
new file mode 100644
index 0000000..2652bfa
--- /dev/null
+++ b/tractvolsx.h
@@ -0,0 +1,217 @@
+/* tractvolsx.h
+
+ Tim Behrens, FMRIB Image Analysis Group
+
+ Copyright (C) 2004 University of Oxford */
+
+/* CCOPYRIGHT */
+
+#ifndef __TRACTVOLSX_H_
+#define __TRACTVOLSX_H_
+
+/////////////////////////////////////////////////////////
+// Class TractVolsx //
+/////////////////////////////////////////////////////////
+
+#include "newimage/newimageall.h"
+#include <iostream>
+#include "stdlib.h"
+#include "probtrackxOptions.h"
+using namespace std;
+using namespace NEWIMAGE;
+using namespace TRACT;
+
+namespace TRACTVOLSX{
+ class Tractvolsx
+ {
+ private:
+ probtrackxOptions& opts;
+ vector<volume4D<float>* > thsamples;
+ vector<volume4D<float>* > phsamples;
+ vector<volume4D<float>* > fsamples;
+ bool init_sample;
+ int fibst;
+ bool usef;
+ public:
+ //constructors::
+ Tractvolsx(const bool& usefin=false):opts(probtrackxOptions::getInstance()),init_sample(true),fibst(1),usef(usefin){}
+ Tractvolsx():opts(probtrackxOptions::getInstance()){}
+ ~Tractvolsx(){
+ for(unsigned int m=0;m<thsamples.size();m++)
+ delete thsamples[m]; //ask flitney, do you just delete the ptr??
+ for(unsigned int m=0;m<phsamples.size();m++)
+ delete phsamples[m];
+ for(unsigned int m=0;m<fsamples.size();m++)
+ delete fsamples[m];
+ }
+
+
+ void reset(const int& fibst_in){
+ init_sample=true;
+ fibst=fibst_in;
+ }
+ //Initialise
+ void initialise(const string& basename){
+
+
+ if(fsl_imageexists(basename+"_thsamples")){
+ volume4D<float> *tmpthptr= new volume4D<float>;
+ volume4D<float> *tmpphptr= new volume4D<float>;
+ volume4D<float> *tmpfptr= new volume4D<float>;
+ cout<<"1"<<endl;
+ read_volume4D(*tmpthptr,basename+"_thsamples");
+ cout<<"2"<<endl;
+ thsamples.push_back(tmpthptr);
+ cout<<"3"<<endl;
+ read_volume4D(*tmpphptr,basename+"_phsamples");
+ cout<<"4"<<endl;
+ phsamples.push_back(tmpphptr);
+ cout<<"5"<<endl;
+ if(usef){
+ read_volume4D(*tmpfptr,basename+"_fsamples");
+ fsamples.push_back(tmpfptr);
+ }
+ cout<<"6"<<endl;
+ }
+ else{
+ int fib=1;
+ bool fib_existed=true;
+ while(fib_existed){
+ if(fsl_imageexists(basename+"_th"+num2str(fib)+"samples")){
+ volume4D<float> *tmpthptr= new volume4D<float>;
+ volume4D<float> *tmpphptr= new volume4D<float>;
+ volume4D<float> *tmpfptr= new volume4D<float>;
+ cout<<fib<<"_1"<<endl;
+ read_volume4D(*tmpthptr,basename+"_th"+num2str(fib)+"samples");
+ thsamples.push_back(tmpthptr);
+ cout<<fib<<"_2"<<endl;
+ read_volume4D(*tmpphptr,basename+"_ph"+num2str(fib)+"samples");
+ phsamples.push_back(tmpphptr);
+ cout<<fib<<"_3"<<endl;
+ read_volume4D(*tmpfptr,basename+"_f"+num2str(fib)+"samples");
+ fsamples.push_back(tmpfptr);
+ fib++;
+ }
+ else{
+ fib_existed=false;
+ }
+ }
+
+ }
+ cout<<"7"<<endl;
+ }
+
+
+ ColumnVector sample(const float& x,const float& y,const float&z,const float& r_x,const float& r_y,const float& r_z){
+
+ ////////Probabilistic interpolation
+ int cx =(int) ceil(x),fx=(int) floor(x);
+ int cy =(int) ceil(y),fy=(int) floor(y);
+ int cz =(int) ceil(z),fz=(int) floor(z);
+
+ //cerr<<x<<" "<<y<<" "<<z<<" "<<cx<<" "<<cy<<" "<<cz<<" "<<fx<<" "<<fy<<" "<<fz<<endl;
+ float pcx,pcy,pcz;
+ if(cx==fx)
+ pcx=1;
+ else
+ pcx=(x-fx)/(cx-fx);
+
+ if(cy==fy)
+ pcy=1;
+ else
+ pcy=(y-fy)/(cy-fy);
+
+ if(cz==fz)
+ pcz=1;
+ else
+ pcz=(z-fz)/(cz-fz);
+
+ ///////new xyz values from probabilistic interpolation
+ int newx,newy,newz;
+ float tmp=rand(); tmp/=RAND_MAX;
+ if(tmp>pcx)
+ newx=fx;
+ else
+ newx=cx;
+
+ tmp=rand(); tmp/=RAND_MAX;
+ if(tmp>pcy)
+ newy=fy;
+ else
+ newy=cy;
+
+ tmp=rand(); tmp/=RAND_MAX;
+ if(tmp>pcz)
+ newz=fz;
+ else
+ newz=cz;
+
+ ColumnVector th_ph_f(3);
+
+
+ float samp=rand(); samp/=RAND_MAX;
+ samp=round(samp*((*thsamples[0]).tsize()-1));
+ float theta=0,phi=0;
+ float dotmax=0,dottmp=0;
+ int fibind=0;
+ if(thsamples.size()>1){//more than 1 fibre
+ if(init_sample){//go for the specified option on the first jump
+ theta=(*thsamples[fibst])(int(newx),int(newy),int(newz),int(samp));
+ phi=(*phsamples[fibst])(int(newx),int(newy),int(newz),int(samp));
+ init_sample=false;
+ }
+ else{
+ for(unsigned int fib=0;fib<thsamples.size();fib++){
+ if((*fsamples[fib])(int(newx),int(newy),int(newz),int(samp))>opts.fibthresh.value()){
+ float phtmp=(*phsamples[fib])(int(newx),int(newy),int(newz),int(samp));
+ float thtmp=(*thsamples[fib])(int(newx),int(newy),int(newz),int(samp));
+ dottmp=fabs(sin(thtmp)*cos(phtmp)*r_x + sin(thtmp)*sin(phtmp)*r_y + cos(thtmp)*r_z);
+ if(dottmp>dotmax){
+ dotmax=dottmp;
+ theta=thtmp;
+ phi=phtmp;
+ fibind=fib;
+ }
+
+ }
+
+ }
+
+ if(dotmax==0){
+ theta=(*thsamples[0])(int(newx),int(newy),int(newz),int(samp));
+ phi=(*phsamples[0])(int(newx),int(newy),int(newz),int(samp));
+ }
+ }
+ }
+ else{
+ theta=(*thsamples[0])(int(newx),int(newy),int(newz),int(samp));
+ phi=(*phsamples[0])(int(newx),int(newy),int(newz),int(samp));
+ }
+
+
+ float f;
+
+ if(usef){
+ f = (*fsamples[fibind])(int(newx),int(newy),int(newz),int(samp));
+ }
+ else
+ f=1;
+
+ th_ph_f(1)=theta;
+ th_ph_f(2)=phi;
+ th_ph_f(3)=f;
+ return th_ph_f;
+ }
+
+ ColumnVector dimensions(){
+ ColumnVector dims(3);
+ dims << (*thsamples[0]).xdim() <<(*thsamples[0]).ydim() << (*thsamples[0]).zdim();
+ return dims;
+ }
+ };
+}
+
+#endif
+
+
+
--
GitLab