diff --git a/ptx_nmasks.cc b/ptx_nmasks.cc index 929bc21e82249d42f57e0497466e266178e77c73..c36f4d0cd18ec49d3745059a48c1e51831fd34c9 100644 --- a/ptx_nmasks.cc +++ b/ptx_nmasks.cc @@ -42,8 +42,7 @@ void nmasks() Seedmanager seedmanager(counter); for(unsigned int i=0;i<seeds.size();i++){ - if(i>0) - stline.pop_waymasks(); + stline.clear_waymasks(); for(unsigned int j=0;j<seeds.size();j++){ if(j!=i) stline.add_waymask(seeds[j]); diff --git a/ptx_seedmask.cc b/ptx_seedmask.cc index c66f0348adc5f9bc38f1b54f29f6771d5e9cf573..d77365d61a01f7efc9328e94b53f2fd6a2d661c0 100644 --- a/ptx_seedmask.cc +++ b/ptx_seedmask.cc @@ -30,19 +30,22 @@ void seedmask() seeds=NEWIMAGE::abs(seeds); seeds.binarise(0,seeds.max()+1,exclusive); + int keeptotal=0; for(int z=0;z<seeds.zsize();z++){ cout <<"sl "<<z<<endl; for(int y=0;y<seeds.ysize();y++){ for(int x=0;x<seeds.xsize();x++){ if(seeds(x,y,z)>0){ cout <<"run"<<endl; - seedmanager.run(x,y,z); + keeptotal += seedmanager.run(x,y,z); } } } } - + + counter.save_total(keeptotal); counter.save(); + cout<<"finished"<<endl; } diff --git a/streamlines.cc b/streamlines.cc index 341aa7e5b87c21fa4eee402cee347155b28d2635..71ea3b74cb6c68ee4b046a2192247ed2986e4e86 100644 --- a/streamlines.cc +++ b/streamlines.cc @@ -214,6 +214,7 @@ namespace TRACT{ } if(rubbish_passed) accept_path=false; + return accept_path; } @@ -548,6 +549,15 @@ namespace TRACT{ } + void Counter::save_total(const int& keeptotal){ + + // save total number of particles that made it through the streamlining + ColumnVector keeptotvec(1); + keeptotvec(1)=keeptotal; + write_ascii_matrix(keeptotvec,logger.appendDir("waytotal")); + + } + void Counter::save(){ if(opts.simpleout.value()){ save_pathdist(); @@ -612,7 +622,8 @@ namespace TRACT{ } - void Seedmanager::run(const float& x,const float& y,const float& z,bool onewayonly, int fibst){ + // this function now returns the total number of pathways that survived a streamlining (SJ) + int Seedmanager::run(const float& x,const float& y,const float& z,bool onewayonly, int fibst){ //onewayonly for mesh things.. cout <<x<<" "<<y<<" "<<z<<endl; if(fibst == -1){ @@ -626,6 +637,7 @@ namespace TRACT{ fibst=1;// fix this for > 2 fibres } + int nlines=0; for(int p=0;p<opts.nparticles.value();p++){ if(opts.verbose.value()>1) logger.setLogFile("particle"+num2str(p)); @@ -637,23 +649,25 @@ namespace TRACT{ forwardflag=true; m_counter.store_path(); m_counter.count_streamline(); + nlines++; } m_stline.reverse(); } if(m_stline.streamline(x,y,z,m_seeddims,fibst)){ backwardflag=true; m_counter.count_streamline(); + //nlines++; //count twice ? } m_counter.clear_streamline(forwardflag,backwardflag); } m_counter.count_seed(); + return nlines; } - } diff --git a/streamlines.h b/streamlines.h index af95f0ab0dbc2dc5e6a3b50f88332d8e384f641b..050937f4d423a2c838dc36f934f378f0cb6d6a8c 100644 --- a/streamlines.h +++ b/streamlines.h @@ -63,6 +63,11 @@ namespace TRACT{ m_own_waymasks.pop_back(); } + void clear_waymasks(){ + // clear all waymasks + for(unsigned int i=0;i<m_waymasks.size();i++) + pop_waymasks(); + } inline const float get_x_seed() const {return m_x_s_init;} inline const float get_y_seed() const {return m_y_s_init;} @@ -165,6 +170,7 @@ namespace TRACT{ void update_maskmatrix(){} //not written yet + void save_total(const int& keeptotal); void save(); void save_pathdist(); void save_pathdist(string add); @@ -197,7 +203,7 @@ namespace TRACT{ m_seeddims.ReSize(3); m_seeddims<<m_seeds.xdim()<<m_seeds.ydim()<<m_seeds.zdim(); } - void run(const float& x,const float& y,const float& z,bool onewayonly=false, int fibst=-1); + int run(const float& x,const float& y,const float& z,bool onewayonly=false, int fibst=-1); }; }