00001
00009 #ifndef GAUSSIANREGIONMODULE_H
00010 #define GAUSSIANREGIONMODULE_H
00011
00012 #include <SegmentationModule/RegionModule.h>
00013 #include <vector>
00014 #include <Schemes.h>
00015
00016 namespace segmentation {
00017
00024 template <unsigned DIM, typename Pixel = float, typename PixelImage = float>
00025 class GaussianRegionModule : public RegionModule<DIM, Pixel, PixelImage>
00026 {
00027 public:
00028
00029 typedef RegionModule<DIM, Pixel, PixelImage> base;
00030 typedef typename base::LevelSet LevelSet;
00031 typedef typename base::Index Index;
00032 typedef typename base::Histogram Histogram;
00033 typedef typename base::ImageType ImageType;
00034
00037
00039 GaussianRegionModule(ImageType * _data,
00040 float _weight = 1.f,
00041 bool _var=false,
00042 bool _adaptive = true
00043 );
00044
00046 virtual ~GaussianRegionModule();
00047
00049
00051
00053 void use_variance(bool _var=true);
00054
00056
00058
00060 float value(const LevelSet & phi, const Index & ind);
00062 float energy(const LevelSet & phi, const Index & ind);
00064 void init(const LevelSet & phi);
00066 void update(const LevelSet & phi);
00068 void save_result(const LevelSet & phi, std::string evodir, int iter);
00069
00070
00072 void update_gaussian_params();
00073
00075
00077
00078 protected:
00079
00080 std::vector<Pixel> m_ll_pos;
00081 std::vector<Pixel> m_ll_neg;
00082 float m_mean_pos;
00083 float m_mean_neg;
00084 float m_var_pos;
00085 float m_var_neg;
00086 bool m_var;
00087
00089
00090 };
00091
00092
00097
00098
00099
00100
00101 template <unsigned DIM, typename Pixel, typename PixelImage>
00102 GaussianRegionModule<DIM, Pixel, PixelImage>::GaussianRegionModule(ImageType * _data, float _weight, bool _var, bool _adaptive)
00103 :RegionModule<DIM, Pixel, PixelImage>(_data, _weight, _adaptive), m_mean_pos(0), m_mean_neg(0), m_var_pos(0), m_var_neg(0), m_var(_var)
00104 { }
00105
00106 template <unsigned DIM, typename Pixel, typename PixelImage>
00107 GaussianRegionModule<DIM, Pixel, PixelImage>::~GaussianRegionModule()
00108 { }
00109
00110
00112
00113
00114
00115 template <unsigned DIM, typename Pixel, typename PixelImage>
00116 void GaussianRegionModule<DIM, Pixel, PixelImage>::use_variance(bool _var)
00117 {
00118 m_var = _var;
00119 }
00120
00121
00123
00124
00125
00126 template <unsigned DIM, typename Pixel, typename PixelImage>
00127 float GaussianRegionModule<DIM, Pixel, PixelImage>::value(const LevelSet & phi, const Index & ind)
00128 {
00129 PixelImage val = (* base::m_data)(ind);
00130 const Pixel ll_pos = m_ll_pos[base::m_hist_pos->get_index(val)];
00131 const Pixel ll_neg = m_ll_neg[base::m_hist_neg->get_index(val)];
00132
00133
00134 return ((m_mean_pos-m_mean_neg < 0)? 1 : -1) * Globals::delta(phi(ind)) * (ll_pos - ll_neg);
00135 }
00136
00137 template <unsigned DIM, typename Pixel, typename PixelImage>
00138 float GaussianRegionModule<DIM, Pixel, PixelImage>::energy(const LevelSet & phi, const Index & ind)
00139 {
00140 PixelImage val = (* base::m_data)(ind);
00141 float ll_pos = m_ll_pos[base::m_hist_pos->get_index(val)];
00142 float ll_neg = m_ll_neg[base::m_hist_neg->get_index(val)];
00143 float h = Globals::heaviside(phi(ind));
00144
00145 return -(1-h)*ll_pos - h*ll_neg;
00146 }
00147
00148 template <unsigned DIM, typename Pixel, typename PixelImage>
00149 void GaussianRegionModule<DIM, Pixel, PixelImage>::init(const LevelSet & phi)
00150 {
00151 base::init(phi);
00152 base::init_histograms(phi);
00153 update_gaussian_params();
00154 }
00155
00156 template <unsigned DIM, typename Pixel, typename PixelImage>
00157 void GaussianRegionModule<DIM, Pixel, PixelImage>::update(const LevelSet & phi)
00158 {
00159 base::update_histograms(phi);
00160 update_gaussian_params();
00161 }
00162
00163 template <unsigned DIM, typename Pixel, typename PixelImage>
00164 void GaussianRegionModule<DIM, Pixel, PixelImage>::update_gaussian_params()
00165 {
00166 m_mean_neg = 0;
00167 m_var_neg = 0;
00168 m_mean_pos = 0;
00169 m_var_pos = 0;
00170
00171
00172 base::hist_neg()->compute_statistics(m_mean_neg, m_var_neg, m_var);
00173 base::hist_pos()->compute_statistics(m_mean_pos, m_var_pos, m_var);
00174
00175
00176 if (m_ll_neg.size() != base::hist_neg()->size()) m_ll_neg.resize(base::hist_neg()->size());
00177 if (m_ll_pos.size() != base::hist_pos()->size()) m_ll_pos.resize(base::hist_pos()->size());
00178
00179 if (m_var) {
00180 for(int i=0; i<m_ll_neg.size(); i++)
00181 m_ll_neg[i] = -0.5f*log(2*M_PI*m_var_neg) - 0.5f*Maths::sqr( i+base::hist_neg()->min() - m_mean_neg)/m_var_neg;
00182 for(int i=0; i<m_ll_pos.size(); i++)
00183 m_ll_pos[i] = -0.5f*log(2*M_PI*m_var_pos) - 0.5f*Maths::sqr( i+base::hist_pos()->min() - m_mean_pos)/m_var_pos;
00184 }
00185 else{
00186 for(int i=0; i<m_ll_neg.size(); i++)
00187 m_ll_neg[i] = - Maths::sqr(i+base::hist_neg()->min() - m_mean_neg);
00188 for(int i=0; i<m_ll_pos.size(); i++)
00189 m_ll_pos[i] = - Maths::sqr(i+base::hist_pos()->min() - m_mean_pos);
00190 }
00191 }
00192
00193 template <unsigned DIM, typename Pixel, typename PixelImage>
00194 void GaussianRegionModule<DIM, Pixel, PixelImage>::save_result(const LevelSet & phi, std::string evodir, int iter)
00195 {
00196
00197
00198 BaseImage<DIM, float> tmp(phi.shape());
00199 std::ostringstream name;
00200
00201 #if 0 // creation image ou regions = couleurs homogene moyenne
00202
00203 for (typename LevelSet::template iterator<domain> i=phi.begin() ; i!=phi.end() ; ++i ){
00204 if(phi(i) > 0) tmp(i) = m_mean_pos ;
00205 else tmp(i) = m_mean_neg ;
00206 }
00207
00208 std::ostringstream name;
00209 std::ofstream file;
00210 std::ofstream file2;
00211
00212 name.str("");
00213 name << evodir << "Gaussian_mean_" << std::setw(6) << std::setfill('0') << iter << ".inr";
00214 file.open(name.str().c_str());
00215 file << format("inrimage") << tmp ;
00216 file.close();
00217 std::cout << " - Saved " << name.str() << " ..." << std::endl;
00218
00219 #endif
00220
00221 #if 0 // image ou regions = couleurs transparente jauen bleu + contour en jaune bleu
00222
00223 typedef Pixels::RGB<unsigned char> RGBPixel;
00224 typedef BaseImage<DIM, RGBPixel> ColorImageType;
00225 ColorImageType tmp2(phi.shape());
00226 for (typename LevelSet::template iterator<domain> i=phi.begin() ; i!=phi.end() ; ++i ){
00227 if(phi(i) > 0) {
00228 tmp2(i) = 0.2*RGBPixel(249,249,0);
00229 const PixelImage _val = (* base::m_data)(i);
00230 PixelImage _min = _val;
00231 _min = ((tmp2(i).green() + _min > 250)? 250 - tmp2(i).green() :_min);
00232 _min = ((tmp2(i).red() + _min > 250)? 250 - tmp2(i).red() :_min);
00233 _min = ((tmp2(i).blue() + _min > 250)? 250 - tmp2(i).blue() :_min);
00234 tmp2(i) += _min;
00235 }
00236 else{
00237 tmp2(i) = 0.2*RGBPixel(40,90,220);
00238 const PixelImage _val = (* base::m_data)(i);
00239 PixelImage _min = _val;
00240 _min = ((tmp2(i).green() + _min > 250)? 250 - tmp2(i).green() :_min);
00241 _min = ((tmp2(i).red() + _min > 250)? 250 - tmp2(i).red() :_min);
00242 _min = ((tmp2(i).blue() + _min > 250)? 250 - tmp2(i).blue() :_min);
00243 tmp2(i) += _min;
00244 }
00245
00246 if(phi.is_in_band(i.position())){
00247 const Index& ind = i.position();
00248 const float _phi = phi(ind);
00249 int inf_coords[DIM];
00250 int sup_coords[DIM];
00251 Schemes<DIM,float>::find_neighbours(phi, ind, inf_coords, sup_coords);
00252
00253 for (int i=0 ; i<DIM ; i++){
00254 Index inf = ind;
00255 Index sup = ind;
00256 inf(i+1) = inf_coords[i];
00257 sup(i+1) = sup_coords[i];
00258 if ((_phi * phi(inf) <= 0) || (_phi * phi(sup) <= 0))
00259 tmp2(ind) = (_phi>0)? RGBPixel(249,249,0): RGBPixel(40,90,220);
00260 }
00261 }
00262 }
00263
00264 name.str("");
00265 name << evodir << "Gaussian_couleur_" << std::setw(6) << std::setfill('0') << iter << ".ppm";
00266 std::ofstream file_pos;
00267 file_pos.open(name.str().c_str());
00268 file_pos << format("rawppm") << tmp2 ;
00269 file_pos.close();
00270 std::cout << " - Saved " << name.str() << " ..." << std::endl;
00271
00272 #endif
00273
00274 #if 0
00275 if (m_var) {
00276
00277 for (typename LevelSet::template iterator<domain> i=phi.begin() ; i!=phi.end() ; ++i ){
00278 if(phi(i) > 0) tmp(i) = m_var_pos;
00279 else tmp(i) = m_var_neg;
00280 }
00281
00282 std::ostringstream name;
00283 std::ofstream file;
00284 std::ofstream file2;
00285
00286 name.str("");
00287 name << evodir << "Gaussian_var_" << std::setw(6) << std::setfill('0') << iter << ".inr";
00288 file.open(name.str().c_str());
00289 file << format("inrimage") << tmp ;
00290 file.close();
00291 std::cout << " - Saved " << name.str() << " ..." << std::endl;
00292
00293
00294
00295 name.str("");
00296 name << evodir << "Gaussian_pos_" << std::setw(6) << std::setfill('0') << iter << ".txt";
00297 file2.open(name.str().c_str());
00298 for(PixelImage i=base::m_hist_pos->min(); i < base::m_hist_pos->max() ; i++)
00299 file2 << i << " " << Globals::gauss_density(m_mean_pos, m_var_pos, i) << std::endl;
00300 file2.close();
00301 std::cout << " - Saved " << name.str() << " ..." << std::endl;
00302
00303 name.str("");
00304 name << evodir << "Gaussian_neg_" << std::setw(6) << std::setfill('0') << iter << ".txt";
00305 file2.open(name.str().c_str());
00306 for(PixelImage i=base::m_hist_neg->min(); i < base::m_hist_neg->max() ; i++)
00307 file2 << i << " " << Globals::gauss_density(m_mean_neg, m_var_neg, i) << std::endl;
00308 file2.close();
00309 std::cout << " - Saved " << name.str() << " ..." << std::endl;
00310 }
00311
00312
00313 #endif
00314 }
00315 }
00316
00317 #endif // GAUSSIANREGIONMODULE_H