From 63c93a94ecbd6bade95e8d87f82f07d5646c4f0a Mon Sep 17 00:00:00 2001 From: Saad Jbabdi <saad@fmrib.ox.ac.uk> Date: Thu, 30 Nov 2006 14:28:20 +0000 Subject: [PATCH] *** empty log message *** --- probtrackxOptions.h | 11 ++++-- streamlines.cc | 94 +++++++++++++++++++++++++++++++++++---------- 2 files changed, 81 insertions(+), 24 deletions(-) diff --git a/probtrackxOptions.h b/probtrackxOptions.h index 36431e5..07724aa 100644 --- a/probtrackxOptions.h +++ b/probtrackxOptions.h @@ -35,13 +35,14 @@ class probtrackxOptions { Option<string> mode; Option<string> targetfile; Option<bool> simpleout; - Option<int> pathdist; + Option<bool> pathdist; Option<bool> s2tout; FmribOption<bool> matrix1out; FmribOption<bool> matrix2out; FmribOption<bool> maskmatrixout; Option<string> outfile; Option<string> rubbishfile; + Option<string> stopfile; Option<string> seeds_to_dti; FmribOption<string> skipmask; Option<string> seedref; @@ -110,9 +111,9 @@ class probtrackxOptions { simpleout(string("--opd"), false, string("output path distribution"), false, no_argument), - pathdist(string("--pd"), 1, + pathdist(string("--pd"), false, string("path distribution, 1:probability,2:distance corrected value"), - false, requires_argument), + false, no_argument), s2tout(string("--os2t"), false, string("output seeds to targets"), false, no_argument), @@ -131,6 +132,9 @@ class probtrackxOptions { rubbishfile(string("--rubbish"), string(""), string("Rubbish file"), false, requires_argument), + stopfile(string("--stop"), string(""), + string("Stop tracking at locations given by this mask file"), + false, requires_argument), seeds_to_dti(string("--xfm"), string(""), string("Transform Matrix taking seed space to DTI space default is to use the identity"),false, requires_argument), skipmask(string("--no_integrity"), string(""), @@ -215,6 +219,7 @@ class probtrackxOptions { options.add(maskmatrixout); options.add(outfile); options.add(rubbishfile); + options.add(stopfile); options.add(seeds_to_dti); options.add(nparticles); options.add(nsteps); diff --git a/streamlines.cc b/streamlines.cc index 330b214..fb42395 100644 --- a/streamlines.cc +++ b/streamlines.cc @@ -37,6 +37,9 @@ namespace TRACT{ if(opts.rubbishfile.value()!=""){ read_volume(m_rubbish,opts.rubbishfile.value()); } + if(opts.stopfile.value()!=""){ + read_volume(m_stop,opts.stopfile.value()); + } vector<string> masknames; if(opts.waypoints.value()!=""){ @@ -90,7 +93,8 @@ namespace TRACT{ x=xst;y=yst;z=zst; m_part.change_xyz(x,y,z); int partlength=0; - //NB - this only goes in one direction!! + bool rubbish_passed=false; + //NB - this only goes in one direction!! for(unsigned int pf=0;pf<m_passed_flags.size();pf++) { m_passed_flags[pf]=false; /// only keep it if this streamline went through all the masks } @@ -139,9 +143,17 @@ namespace TRACT{ if(opts.rubbishfile.value()!=""){ - if(m_rubbish(x_s,y_s,z_s)>0) break; + if(m_rubbish(x_s,y_s,z_s)!=0){ + rubbish_passed=true; + break; + } } - + if(opts.stopfile.value()!=""){ + if(m_stop(x_s,y_s,z_s)!=0){ + m_path.pop_back(); + break; + } + } if(opts.skipmask.value() == ""){ @@ -189,6 +201,8 @@ namespace TRACT{ if(!m_passed_flags[i]) accept_path=false; } + if(rubbish_passed) + accept_path=false; return accept_path; } @@ -363,7 +377,7 @@ namespace TRACT{ void Counter::update_pathdist(){ const vector<ColumnVector>& path=m_stline.get_path_ref(); - if(opts.pathdist.value()==1){ + if(!opts.pathdist.value()){ for(unsigned int i=0;i<path.size();i++){ int x_s=int(round(float(path[i](1)))),y_s=int(round(float(path[i](2)))),z_s=int(round(float(path[i](3)))); if(m_beenhere(x_s,y_s,z_s)==0){ @@ -407,13 +421,34 @@ namespace TRACT{ int xseedvox=int(round(m_stline.get_x_seed())); int yseedvox=int(round(m_stline.get_y_seed())); int zseedvox=int(round(m_stline.get_z_seed())); - for(unsigned int i=0;i<path.size();i++){ - int x_s=int(round(float(path[i](1)))),y_s=int(round(float(path[i](2)))),z_s=int(round(float(path[i](3)))); - for(unsigned int m=0;m<m_targetmasknames.size();m++){ - if(m_targetmasks[m](x_s,y_s,z_s)>0 && m_targflags[m]==0){ - m_seedcounts[m](xseedvox,yseedvox,zseedvox)=m_seedcounts[m](xseedvox,yseedvox,zseedvox)+1; - m_targflags[m]=1; - //m_particle_numbers[m].push_back(particle_number); + + if(!opts.pathdist.value()){ + for(unsigned int i=0;i<path.size();i++){ + int x_s=int(round(float(path[i](1)))),y_s=int(round(float(path[i](2)))),z_s=int(round(float(path[i](3)))); + for(unsigned int m=0;m<m_targetmasknames.size();m++){ + if(m_targetmasks[m](x_s,y_s,z_s)>0 && m_targflags[m]==0){ + m_seedcounts[m](xseedvox,yseedvox,zseedvox)=m_seedcounts[m](xseedvox,yseedvox,zseedvox)+1; + m_targflags[m]=1; + //m_particle_numbers[m].push_back(particle_number); + } + } + } + } + else{ + int d=0; + int x_s=0,y_s=0,z_s=0,px,py,pz; + for(unsigned int i=0;i<path.size();i++){ + px=(int)floor(x_s),py=(int)floor(y_s),pz=(int)floor(z_s); + x_s=int(round(float(path[i](1))));y_s=int(round(float(path[i](2))));z_s=int(round(float(path[i](3)))); + if(i==0 | (floor(x_s)!=px | floor(y_s)!=py | floor(z_s)!=pz)){ + d++; + for(unsigned int m=0;m<m_targetmasknames.size();m++){ + if(m_targetmasks[m](x_s,y_s,z_s)>0 && m_targflags[m]==0){ + m_seedcounts[m](xseedvox,yseedvox,zseedvox)+=d; + m_targflags[m]=1; + //m_particle_numbers[m].push_back(particle_number); + } + } } } } @@ -446,18 +481,35 @@ namespace TRACT{ void Counter::update_matrix2_row(){ //run this one every streamline - not every voxel.. const vector<ColumnVector>& path=m_stline.get_path_ref(); - for(unsigned int i=0;i<path.size();i++){ - ColumnVector xyz_seeds=path[i]; - ColumnVector xyz_lr=vox_to_vox(xyz_seeds,m_seedsdim,m_lrdim,m_I); - int x_lr=int(round(float(xyz_lr(1)))),y_lr=int(round(float(xyz_lr(2)))),z_lr=int(round(float(xyz_lr(3)))); - int Concol2=m_lookup2(x_lr,y_lr,z_lr,0); - if(Concol2!=0){ - if(m_beenhere2(x_lr,y_lr,z_lr)==0){ - m_ConMat2(m_Conrow2,Concol2,0)+=1; - m_beenhere2(x_lr,y_lr,z_lr)=1; + + if(!opts.pathdist.value()) + for(unsigned int i=0;i<path.size();i++){ + ColumnVector xyz_seeds=path[i]; + ColumnVector xyz_lr=vox_to_vox(xyz_seeds,m_seedsdim,m_lrdim,m_I); + int x_lr=int(round(float(xyz_lr(1)))),y_lr=int(round(float(xyz_lr(2)))),z_lr=int(round(float(xyz_lr(3)))); + int Concol2=m_lookup2(x_lr,y_lr,z_lr,0); + if(Concol2!=0){ + if(m_beenhere2(x_lr,y_lr,z_lr)==0){ + m_ConMat2(m_Conrow2,Concol2,0)+=1; + m_beenhere2(x_lr,y_lr,z_lr)=1; + } + } + + } + else{ + int d=1; + for(unsigned int i=0;i<path.size();i++){ + ColumnVector xyz_seeds=path[i]; + ColumnVector xyz_lr=vox_to_vox(xyz_seeds,m_seedsdim,m_lrdim,m_I); + int x_lr=int(round(float(xyz_lr(1)))),y_lr=int(round(float(xyz_lr(2)))),z_lr=int(round(float(xyz_lr(3)))); + int Concol2=m_lookup2(x_lr,y_lr,z_lr,0); + if(Concol2!=0){ + if(m_beenhere2(x_lr,y_lr,z_lr)==0){ + m_ConMat2(m_Conrow2,Concol2,0)+=d;d++; + m_beenhere2(x_lr,y_lr,z_lr)=1; + } } } - } } -- GitLab