iRoCS Toolbox  1.1.0
iRoCS.hh
Go to the documentation of this file.
1 /**************************************************************************
2  *
3  * Copyright (C) 2015 Thorsten Falk
4  *
5  * Image Analysis Lab, University of Freiburg, Germany
6  *
7  * This program is free software; you can redistribute it and/or modify
8  * it under the terms of the GNU General Public License as published by
9  * the Free Software Foundation; either version 3 of the License, or
10  * (at your option) any later version.
11  *
12  * This program is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  * GNU General Public License for more details.
16  *
17  * You should have received a copy of the GNU General Public License
18  * along with this program; if not, write to the Free Software Foundation,
19  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
20  *
21  **************************************************************************/
22 
23 /*======================================================================*/
32 /*======================================================================*/
33 
34 #ifndef ATBIROCS_HH
35 #define ATBIROCS_HH
36 
37 #ifdef HAVE_CONFIG_H
38 #include <config.hh>
39 #endif
40 
41 #include <vector>
42 #include <blitz/array.h>
43 
46 
48 #include "ATBNucleus.hh"
49 #include "Array.hh"
50 #include "Interpolator.hh"
51 
52 namespace atb
53 {
54 
55 /*======================================================================*/
68 /*======================================================================*/
69  class IRoCS
70  {
71 
72  public:
73 
74 /*======================================================================*/
83 /*======================================================================*/
84  IRoCS(iRoCS::ProgressReporter *progressReporter = NULL);
85 
86 /*======================================================================*/
90 /*======================================================================*/
91  ~IRoCS();
92 
93 /*======================================================================*/
101 /*======================================================================*/
102  void setProgressReporter(iRoCS::ProgressReporter *progressReporter);
103 
104 /*======================================================================*/
146 /*======================================================================*/
147  void fit(blitz::TinyVector<double,3> const &qcPositionUm,
148  std::list<Nucleus> const &nuclei,
149  double kappa = 1.0, double lambda = 0.0, double mu = 0.0,
150  int nIter = 1000, double tau = 0.1, double searchRadiusUm = 0.0);
151 
152 /*======================================================================*/
194 /*======================================================================*/
195  void fit(blitz::TinyVector<double,3> const &qcPositionUm,
196  std::vector<Nucleus> const &nuclei,
197  double kappa = 1.0, double lambda = 0.0, double mu = 0.0,
198  int nIter = 1000, double tau = 0.1, double searchRadiusUm = 0.0);
199 
200 /*======================================================================*/
242 /*======================================================================*/
243  void fit(blitz::TinyVector<double,3> const &qcPositionUm,
244  Array<int,3> const &segmentation,
245  double kappa = 1.0, double lambda = 0.0, double mu = 0.0,
246  int nIter = 1000, double tau = 0.1, double searchRadiusUm = 0.0);
247 
248 /*======================================================================*/
289 /*======================================================================*/
290  void fit(blitz::TinyVector<double,3> const &qcPositionUm,
291  std::vector< blitz::TinyVector<double,3> > &positionsUm,
292  double kappa = 1.0, double lambda = 0.0, double mu = 0.0,
293  int nIter = 1000, double tau = 0.1, double searchRadiusUm = 0.0);
294 
295 /*======================================================================*/
305 /*======================================================================*/
306  blitz::TinyVector<double,3> getCoordinates(
307  blitz::TinyVector<double,3> const &pos) const;
308 
309 /*======================================================================*/
320 /*======================================================================*/
321  double getDistanceToSurface(blitz::TinyVector<double,3> const &pos) const;
322 
323 /*======================================================================*/
334 /*======================================================================*/
335  blitz::TinyVector<double,3> getEuclideanPositionUm(
336  blitz::TinyVector<double,3> const &coordinates) const;
337 
338 /*======================================================================*/
352 /*======================================================================*/
353  template<typename DataT>
355  blitz::Array<DataT,3> const &data,
356  blitz::TinyVector<double,3> const &originalElementSizeUm,
357  blitz::Array<DataT,3> &straightened,
358  blitz::TinyVector<double,3> const &straightenedElementSizeUm,
359  blitz::TinyVector<double,3> const &originUm, double phiOffset) const;
360 
361  blitz::TinyVector<double,3> getAxisPosition(double u) const;
362  double thickness(double u) const;
363 
364  CoupledBSplineModel<3> const &ccm() const;
366  BSpline<double> const &thicknessSpline() const;
367 
368  double uQC() const;
369  blitz::TinyMatrix<double,4,4> const &normalizationTransformation() const;
370  blitz::TinyMatrix<double,4,4> const &inverseNormalizationTransformation()
371  const;
372 
373  blitz::Array<blitz::TinyVector<double,3>,1> &vertices();
374  blitz::Array<blitz::TinyVector<ptrdiff_t,3>,1> &indices();
375  blitz::Array<blitz::TinyVector<double,3>,1> &normals();
376 
377  void save(std::string const &fileName, std::string const &groupName,
378  bool throwErrors = false) const;
379  void load(std::string const &fileName, std::string const &groupName,
380  bool throwErrors = false);
381 
382  void save(BlitzH5File &outFile, std::string const &groupName) const;
383  void load(BlitzH5File const &inFile, std::string const &groupName);
384 
385  private:
386 
387  void extractBoundary();
388  void computeTransformation();
389 
390  void fitEpidermalQuadricRANSAC();
391 
392  iRoCS::ProgressReporter* p_progress;
393 
394  CoupledBSplineModel<3> *p_ccm;
395  double _kappa, _lambda, _mu, _tau, _searchRadiusUm;
396  int _nIter;
397 
398  blitz::TinyMatrix<double,4,4> _trafo, _trafoInv;
399  blitz::TinyVector<double,3> _qcPos;
400  double _uQC;
401 
402  ptrdiff_t _nLatitudes, _nLongitudes;
403  blitz::Array<blitz::TinyVector<double,3>,1> _vertices;
404  blitz::Array<blitz::TinyVector<ptrdiff_t,3>,1> _indices;
405  blitz::Array<blitz::TinyVector<double,3>,1> _normals;
406 
407  double* p_curveLengthCache;
408 
409  };
410 
411 }
412 
413 #include "iRoCS.icc"
414 
415 #endif
void save(std::string const &fileName, std::string const &groupName, bool throwErrors=false) const
Array class derived from blitz++ Arrays for handling microscopic datasets with associated element siz...
The BSpline class provides functions for fitting B-Splines to point clouds and evaluating them at arb...
Definition: ATBSpline.hh:68
The IRoCS class provides means to attach iRoCS to different kinds of processed root images...
Definition: iRoCS.hh:69
blitz::TinyVector< double, 3 > getCoordinates(blitz::TinyVector< double, 3 > const &pos) const
Get the iRoCS position for the given image position in micrometers.
double thickness(double u) const
BSpline< double > const & thicknessSpline() const
Lightweight alternative to libBlitzHDF5 providing its basic functionality.
blitz::Array< blitz::TinyVector< ptrdiff_t, 3 >, 1 > & indices()
BSpline< blitz::TinyVector< double, 3 > > const & axisSpline() const
double getDistanceToSurface(blitz::TinyVector< double, 3 > const &pos) const
Get the signed euclidean distance of the given point to the tubular surface.
Nucleus class containing cell nucelus specific parameters.
void computeStraightenedView(blitz::Array< DataT, 3 > const &data, blitz::TinyVector< double, 3 > const &originalElementSizeUm, blitz::Array< DataT, 3 > &straightened, blitz::TinyVector< double, 3 > const &straightenedElementSizeUm, blitz::TinyVector< double, 3 > const &originUm, double phiOffset) const
Compute a straightened view of the given dataset according to this iRoCS instance.
IRoCS(iRoCS::ProgressReporter *progressReporter=NULL)
Default constructor.
void setProgressReporter(iRoCS::ProgressReporter *progressReporter)
Attach a iRoCS::ProgressReporter through which progress will be output.
Classes and functions for sub-pixel Array access with different interpolation strategies.
~IRoCS()
Destructor.
blitz::Array< blitz::TinyVector< double, 3 >, 1 > & normals()
blitz::TinyMatrix< double, 4, 4 > const & normalizationTransformation() const
blitz::TinyVector< double, 3 > getAxisPosition(double u) const
void load(std::string const &fileName, std::string const &groupName, bool throwErrors=false)
blitz::TinyVector< double, 3 > getEuclideanPositionUm(blitz::TinyVector< double, 3 > const &coordinates) const
Get the euclidean position in micrometers for the given iRoCS coordinate vector.
blitz::Array< blitz::TinyVector< double, 3 >, 1 > & vertices()
double uQC() const
void fit(blitz::TinyVector< double, 3 > const &qcPositionUm, std::list< Nucleus > const &nuclei, double kappa=1.0, double lambda=0.0, double mu=0.0, int nIter=1000, double tau=0.1, double searchRadiusUm=0.0)
Fit iRoCS to a given list of nuclei.
blitz::TinyMatrix< double, 4, 4 > const & inverseNormalizationTransformation() const
CoupledBSplineModel< 3 > const & ccm() const