A.I.M.S algorithms


roistat.h
Go to the documentation of this file.
1 /* This software and supporting documentation are distributed by
2  * Institut Federatif de Recherche 49
3  * CEA/NeuroSpin, Batiment 145,
4  * 91191 Gif-sur-Yvette cedex
5  * France
6  *
7  * This software is governed by the CeCILL-B license under
8  * French law and abiding by the rules of distribution of free software.
9  * You can use, modify and/or redistribute the software under the
10  * terms of the CeCILL-B license as circulated by CEA, CNRS
11  * and INRIA at the following URL "http://www.cecill.info".
12  *
13  * As a counterpart to the access to the source code and rights to copy,
14  * modify and redistribute granted by the license, users are provided only
15  * with a limited warranty and the software's author, the holder of the
16  * economic rights, and the successive licensors have only limited
17  * liability.
18  *
19  * In this respect, the user's attention is drawn to the risks associated
20  * with loading, using, modifying and/or developing or reproducing the
21  * software by the user in light of its specific status of free software,
22  * that may mean that it is complicated to manipulate, and that also
23  * therefore means that it is reserved for developers and experienced
24  * professionals having in-depth computer knowledge. Users are therefore
25  * encouraged to load and test the software's suitability as regards their
26  * requirements in conditions enabling the security of their systems and/or
27  * data to be ensured and, more generally, to use and operate it in the
28  * same conditions as regards security.
29  *
30  * The fact that you are presently reading this means that you have had
31  * knowledge of the CeCILL-B license and that you accept its terms.
32  */
33 
34 
35 #ifndef AIMS_ROI_ROISTAT_H
36 #define AIMS_ROI_ROISTAT_H
37 
39 #include <aims/roi/roi.h>
40 #include <aims/roi/roistat.h>
41 #include <aims/roi/roiselector.h>
42 #include <aims/roi/roigtm.h>
43 #include <aims/data/pheader.h>
45 #include <aims/math/mathelem.h>
46 #include <graph/graph/graph.h>
47 #include <float.h>
48 
49 namespace aims
50 {
51  class RoiSelector;
52  class RoiGtm;
53 
54  class AIMSALGO_API RoiStats : public Graph
55  {
56  public:
57 
58  RoiStats( Resampler< float >* interpolator , RoiSelector* roiSel );
59  virtual ~RoiStats( );
60 
61  template <typename T>
62  void populate( AimsData<T> &image, AimsRoi &roi,
63  const Motion &motion );
64  std::vector<float> tacByStructName(std::string s);
65 
66  void applyGtm( RoiGtm& rgtm );
67  void streamout(char *out=NULL);
68  bool isGtmCorrected() { return _isGtmCorrected; }
69  void setGtmCorrected() { _isGtmCorrected = true; }
70 
71  protected:
75 
76  };
77 
78 
79  template <typename T>
80  void RoiStats::populate( AimsData<T> &image, AimsRoi &roi,
81  const Motion &motion )
82  {
83  std::set< Vertex*> currentRoiStatVert; //tmp var for current set vert of this
84  Vertex* v; //tmp var for vertice of the object "this"
85  std::vector< int > roi_st, roi_dt;
86 
87  // Tatoo the roi object AND the RoiStat object with the image name
88  roi.setProperty( "last_image_stated", image);
89  setProperty("last_image_stated", image);
90 
91  // get the maximum bouding box
92  std::vector< int > bbmaxTO;
93  std::vector< float > vsTO;
94  roi.getProperty( "boundingbox_max", bbmaxTO);
95  roi.getProperty( "voxel_size", vsTO );
96 
97  // Gestion de start et duration compte tenu d'appel monoframe ou multi frame
98  std::vector< int > i_st, i_dt;
99  PythonHeader *hd = dynamic_cast<PythonHeader *>( image.header() );
100  if (! hd) {
101  for(int i= 0; i < image.dimT(); i++) {
102  i_st.push_back(i);
103  i_dt.push_back(1);
104  }
105  } else {
106  if ( hd->hasProperty("start_time") ) {
107  hd->getProperty("start_time", i_st);
108  } else {
109  ASSERT(image.dimT() == 1); // dynamic series with no start/dur time
110  i_st.push_back(0); // Force val for "static" series
111  }
112  if ( hd->hasProperty("duration_time") ) {
113  hd->getProperty("duration_time", i_dt);
114  } else {
115  ASSERT(image.dimT() == 1); // dynamic series with no start/dur time
116  i_dt.push_back(1); // Force val for "static" series
117  }
118  }
119 
120  // Start gathering data.
121  AimsData< float > tmpin( image.dimX(), image.dimY(), image.dimZ() );
122  tmpin.setSizeXYZT( image.sizeX(), image.sizeY(), image.sizeZ(), 1.0 );
123  if(image.header() != 0)
124  tmpin.setHeader(image.header()->cloneHeader()) ;
125  const std::list<BaseTree*> & metaRoiList = _roiSel->children();
128 
129  // Level 1 : Start loop on time dimension
130  for (int t = 0; t < image.dimT(); ++t) // populate may be invoked in byframe
131  { // or allframe mode
132  // 1-Init des fields start et duration time;
133  // 2-Re Copy of one frame in tem image
134  // 3-Set in order the interpolation process
135  getProperty( "start_time", roi_st );
136  roi_st.push_back( i_st[ t ] );
137  setProperty("start_time",roi_st);
138  getProperty( "duration_time", roi_dt );
139  roi_dt.push_back( i_dt[ t ] );
140  setProperty("duration_time",roi_dt);
141 
142  int x,y,z;
143  ForEach3d( tmpin, x, y, z )
144  tmpin(x, y, z) = image(x, y, z, t);
145 
146  _interpolator->setRef( tmpin );
147  AimsData< float > tmpout = _interpolator->doit( motion,
148  bbmaxTO[0]+1, bbmaxTO[1]+1, bbmaxTO[2]+1,
149  Point3df( vsTO[0],vsTO[1],vsTO[2] ) );
150 
151  // Level 2 : Loop on the different MetaRoi of selectionSet
152  for (BaseTree::const_iterator metaRoi = metaRoiList.begin();
153  metaRoi != metaRoiList.end();
154  metaRoi++)
155  {
156  std::map<int,carto::rc_ptr<BucketMap<Void> > > corresBuckets;
157  std::map<int,carto::rc_ptr<BucketMap<Void> > >::iterator
158  ib, eb = corresBuckets.end();
159  Tree* tmpTree = dynamic_cast< Tree *>( *metaRoi ); // Test
160  std::string surname;(*tmpTree).getProperty("surname", surname);
161 
162  // Construct the corresBucket
163  const std::list<BaseTree*> & metaRoiChild = (*metaRoi)->children();
164  for( BaseTree::const_iterator metaRoiContent = metaRoiChild.begin();
165  metaRoiContent != metaRoiChild.end();
166  ++metaRoiContent )
167  {
168  tmpTree = dynamic_cast< Tree *>( *metaRoiContent ); // Test
169  std::string tmp;(*tmpTree).getProperty("nomenName", tmp);
170  //std::cout << "Debug nomenName>> " << tmp << std::endl;
171 
172  std::set<Vertex*> roisv = roi.VertByNameAndDescendantName( tmp );
173  for( std::set< Vertex* >::const_iterator j = roisv.begin();
174  j != roisv.end(); ++j )
175  {
176  int l; (*j)->getProperty("roi_label", l);
177  if( (*j)->getProperty( "aims_roi", bck ) )
178  corresBuckets[ l ] = bck;
179  else
180  std::cout << "Bucket not found. Label = " << l << std::endl ;
181  }
182  }
183 
184  if (corresBuckets.size() == 0) continue;
185 
186  // Lets process the current corresBuckets
187  // and construct the vertices of RoiStat.
188  currentRoiStatVert = getVerticesWith("surname", surname);
189 
190  ASSERT (currentRoiStatVert.size() < 2 ); // only 0 first loop
191  //or 1 -other loops
192  if (currentRoiStatVert.size() == 0)
193  {
194  v = addVertex();
195  v->setProperty("surname", surname);
196  }
197  else v = *( currentRoiStatVert.begin() );
198 
199  // Level 3: loop on the corresBuckets in the arg roi and accumulate
200  float mean = 0.0;
201  float std = 0.0;
202  float min = FLT_MAX;
203  float max = FLT_MIN;
204  int pn2 = 0 ;
205  for( ib=corresBuckets.begin(); ib!=eb; ++ib )
206  {
207  bck = ib->second;
208  // 0- Init
209  // 1-Retrieve one bucket : (*v)
210  // 2-Cumul stat
211  if( !bck->empty() )
212  for( ibk=bck->begin()->second.begin(),
213  ebk=bck->begin()->second.end(); ibk!=ebk; ++ibk )
214  {
215  float pixel = tmpout( ibk->first );
216  mean += pixel;
217  std += square( pixel);
218  min = (min > pixel) ? pixel : min;
219  max = (max < pixel) ? pixel : max;
220  ++pn2 ;
221  }
222  } // end of Level 3 loop
223 
224  //Update the current vertice in the RoiStat object
225  mean /= (float) pn2;
226  std=(float)sqrt( (double)(std-square(mean)*pn2) / (float)(pn2-1) );
227 
228  v->setProperty("voxel_number", pn2 );
229  std::vector< float > tac, sac, mac, Mac;
230  v->getProperty( "mean_ac", tac);
231  tac.push_back( mean );
232  v->setProperty("mean_ac", tac);
233  v->getProperty( "std_ac", sac);
234  sac.push_back( std );
235  v->setProperty( "std_ac", sac);
236  v->getProperty( "min_ac", mac);
237  mac.push_back( min );
238  v->setProperty( "min_ac", mac);
239  v->getProperty( "max_ac", Mac);
240  Mac.push_back( max );
241  v->setProperty( "max_ac", Mac);
242  } //end of level 2 (loop on different set in the selectionSet)
243  } // end level 1 (loop on time)
244 
245  setProperty("voxel_volume",vsTO[0]*vsTO[1]*vsTO[2]);
246  }
247 }
248 
249 #endif
float min(float x, float y)
Definition: thickness.h:105
float max(float x, float y)
Definition: thickness.h:97
Vertex * addVertex(const std::string &s="")
int dimX() const
float sizeY() const
virtual bool hasProperty(const std::string &) const
#define ForEach3d(thing, x, y, z)
int dimT() const
void setGtmCorrected()
Definition: roistat.h:69
std::set< Vertex * > getVerticesWith(const std::string &s) const
void setRef(const AimsData< T > &ref)
set the input data to be resampled
Definition: resampler.h:125
std::list< BaseTree * >::const_iterator const_iterator
virtual Header * cloneHeader(bool keepUuid=false) const =0
std::set< Vertex * > VertByNameAndDescendantName(std::string &name)
Resampler< float > * _interpolator
Definition: roistat.h:72
VSet::iterator iterator
const std::list< BaseTree * > & children() const
virtual void doit(const aims::AffineTransformation3d &transform, AimsData< T > &output_data)
Resample the reference input data (set via setRef()) into an existing output data.
Definition: resampler_d.h:169
RoiSelector * _roiSel
Definition: roistat.h:73
float sizeZ() const
bool isGtmCorrected()
Definition: roistat.h:68
virtual bool getProperty(const std::string &, Object &) const
float sizeX() const
int dimZ() const
void setSizeXYZT(float sizex=1.0f, float sizey=1.0f, float sizez=1.0f, float sizet=1.0f)
int dimY() const
void populate(AimsData< T > &image, AimsRoi &roi, const Motion &motion)
Definition: roistat.h:80
#define AIMSALGO_API
#define ASSERT(EX)
const aims::Header * header() const
T square(const T &val)
bool _isGtmCorrected
Definition: roistat.h:74