00001 #ifndef MOTIONPROFILE_H
00002 #define MOTIONPROFILE_H 1
00003
00004 #include "ncutSetting.h"
00005 #include "matrix.h"
00006 #include "ncutImage.h"
00007 #include "ncutException.h"
00008
00009
00010 #include <math.h>
00011 #include <cstdlib>
00012 #include <iostream>
00013
00014
00015
00016
00017 #define MAX_SEG_ID 2000000000
00018
00019 namespace ncut
00020 {
00021
00022 const double PI = 3.14159265;
00023
00024 struct pxlInfo_t {
00025 double hueVecX;
00026 double hueVecY;
00027 double lit;
00028 unsigned int elm;
00029 };
00030
00044 template <class IMG_ELM_TYPE = unsigned char,
00045 class MSK_ELM_TYPE = signed long>
00046 class Profile
00047 {
00048 public:
00049
00050
00051
00056 inline const Matrix* similarity() const
00057 {return &similarity_;};
00062 inline const std::map<MSK_ELM_TYPE, unsigned int>* indices() const
00063 {return &indices_;};
00068 inline const std::vector<MSK_ELM_TYPE>* identifiers() const
00069 {return &identifiers_;};
00075 unsigned int idx(const MSK_ELM_TYPE& id) const;
00081 inline const MSK_ELM_TYPE& id(unsigned int idx) const
00082 {return identifiers_[idx];};
00083
00084
00085
00092 Profile(const Setting* setting);
00093
00099 Profile(const Profile& clone);
00100
00112 Profile(const Matrix& similarity,
00113 const std::map<MSK_ELM_TYPE, unsigned int>& indices,
00114 const std::vector<MSK_ELM_TYPE>& identifiers,
00115 const Setting* setting);
00116
00120 virtual ~Profile();
00121
00128 virtual Profile& operator=(const Profile& clone);
00129
00130 protected:
00131
00132 Matrix similarity_;
00136 std::map<MSK_ELM_TYPE, unsigned int> indices_;
00140 std::vector<MSK_ELM_TYPE> identifiers_;
00143 const Setting* setting_;
00144 };
00145
00158 template <class IMG_ELM_TYPE = unsigned char,
00159 class MSK_ELM_TYPE = signed long>
00160 class ElementProfile : public Profile<IMG_ELM_TYPE, MSK_ELM_TYPE>
00161 {
00162 public:
00163
00164
00165 const double normFac;
00168
00169
00175 inline const Image<IMG_ELM_TYPE>* frame(unsigned int num) const
00176 {return frameSeq_[num];};
00182 inline const Image<MSK_ELM_TYPE>* mask(unsigned int num) const
00183 {return maskSeq_[num];};
00190 inline const MSK_ELM_TYPE& mask(unsigned int num, unsigned int pxl) const
00191 {return maskSeq_[num]->val(pxl);};
00192
00199 inline const IMG_ELM_TYPE& frame(unsigned int num, unsigned int pxlelm) const
00200 {return frameSeq_[num]->val(pxlelm);};
00201
00206 inline unsigned int width() const
00207 {return width_;};
00212 inline unsigned int height() const
00213 {return height_;};
00218 inline unsigned int nChl() const
00219 {return nChl_;};
00225 inline unsigned int nFrm() const
00226 {return frameSeq_.length();};
00231 inline unsigned int wSize() const
00232 {return wSize_;};
00237 inline unsigned int nElm() const
00238 {return nElm_;};
00245 inline const std::vector<std::vector<unsigned int> >* pixels() const
00246 {return &pxls_;};
00252 inline const std::vector<unsigned int>* pixels(unsigned int idx) const
00253 {return &pxls_[idx];};
00259 unsigned int frameOf(unsigned int idx) const;
00264 inline const Matrix* neighbors() const
00265 {return &neighbors_;};
00271 inline double hueVecX(unsigned int elm) const
00272 {return hueVecX_[elm];};
00278 inline double hueVecY(unsigned int elm) const
00279 {return hueVecY_[elm];};
00285 inline double lightness(unsigned int elm) const
00286 {return lightness_[elm];};
00292 inline double saturation(unsigned int elm) const
00293 {return saturation_[elm];};
00300 inline double intensity(unsigned int elm, unsigned int chl) const
00301 {return intensity_[elm][chl];};
00307 inline double locationX(unsigned int elm) const
00308 {return locationX_[elm];};
00314 inline double locationY(unsigned int elm) const
00315 {return locationY_[elm];};
00323 inline double locationZ(unsigned int elm) const
00324 {return locationZ_[elm];};
00330 inline unsigned int nPixels(unsigned int elm) const
00331 {return nPixels_[elm];};
00337 inline unsigned int minIdx(unsigned int frm) const
00338 {return (frm == 0) ? 0 : maxElmIdx_[frm-1]+1;};
00344 inline unsigned int maxIdx(unsigned int frm) const
00345 {return maxElmIdx_[frm];};
00346
00347
00348
00359 int push(const Image<IMG_ELM_TYPE>* frame, const Image<MSK_ELM_TYPE>* mask);
00360
00361
00362
00372 ElementProfile(unsigned int wSize, const Setting* setting);
00373
00379 ElementProfile(const ElementProfile& clone);
00380
00384 virtual ~ElementProfile();
00385
00391 virtual ElementProfile& operator=(const ElementProfile& clone);
00392
00393 protected:
00394
00395 Sequence<IMG_ELM_TYPE> frameSeq_;
00397 Sequence<MSK_ELM_TYPE> maskSeq_;
00400 unsigned int wSize_;
00402 unsigned int width_;
00404 unsigned int height_;
00406 unsigned int nChl_;
00408 unsigned int nElm_;
00411 std::vector<std::vector<unsigned int> > pxls_;
00414 std::vector<unsigned int> maxElmIdx_;
00417 Matrix neighbors_;
00420
00421 double* hue_;
00422 double* hueVariation_;
00423 double* hueVecX_;
00424 double* hueVecY_;
00425 double* lightness_;
00426 double* litVariation_;
00427 double* saturation_;
00428 double* satVariation_;
00429 double** intensity_;
00430 double* locationX_;
00431 double* locationY_;
00432 double* locationZ_;
00433 unsigned int* nPixels_;
00435 private:
00436
00443 void readaptBase();
00452 void chartElements(unsigned int frmNo);
00456 void initInfoArrays();
00460 void initInfo();
00467 void readaptInfo();
00477 void getInfo(unsigned int from, unsigned int to);
00487 void calcSim(unsigned int from, unsigned int to);
00488
00492 void cleanup();
00493
00494 unsigned int nElm_old_;
00496 unsigned int nElm_front_;
00499 };
00500
00512 template <class IMG_ELM_TYPE = unsigned char,
00513 class MSK_ELM_TYPE = signed long>
00514 class SegmentProfile : public Profile<IMG_ELM_TYPE, MSK_ELM_TYPE>
00515 {
00516 public:
00517
00518
00519
00526 inline const std::vector<std::vector<MSK_ELM_TYPE> >* elements() const
00527 {return &elements_;};
00533 inline const std::vector<MSK_ELM_TYPE>* elements(unsigned int idx) const
00534 {return &elements_[idx];};
00539 inline const ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>* elmProfile() const
00540 {return elmProfile_;};
00546 const MSK_ELM_TYPE& mask(const MSK_ELM_TYPE& elmId) const;
00552 inline unsigned int idxMask(unsigned int elmIdx) const
00553 {return idxMask_[elmIdx];};
00559 unsigned int age(const MSK_ELM_TYPE& segId) const;
00565 bool inactive(const MSK_ELM_TYPE& segId) const;
00571 inline const std::map<MSK_ELM_TYPE, bool>* inactive() const
00572 {return &inactive_;};
00577 inline const std::map<MSK_ELM_TYPE, unsigned int>* age() const
00578 {return &age_;};
00584 const MSK_ELM_TYPE& simSeg(const MSK_ELM_TYPE& segId) const;
00589 inline unsigned int nElm() const
00590 {return (elmProfile_==NULL)?0:elmProfile_->nElm();};
00595 inline unsigned int nSeg() const
00596 {return nSeg_;};
00601 inline unsigned int nChl() const
00602 {return nChl_;};
00608 inline double locationX(unsigned int seg) const
00609 {return locationX_[seg];};
00615 inline double locationY(unsigned int seg) const
00616 {return locationY_[seg];};
00622 inline unsigned int nPixels(unsigned int seg) const
00623 {return nPixels_[seg];};
00624
00625
00626
00639 double
00640 dif(unsigned int idx1, unsigned int idx2,
00641 const std::vector<std::vector<MSK_ELM_TYPE> >* elements = NULL) const;
00642
00653 int push(const ElementProfile<IMG_ELM_TYPE,MSK_ELM_TYPE>* elmProfile,
00654 const std::vector<std::vector<MSK_ELM_TYPE> >& segments);
00655
00663 int calcSim();
00664
00665
00666
00674 SegmentProfile(const Setting* setting);
00675
00681 SegmentProfile(const SegmentProfile& clone);
00682
00686 virtual ~SegmentProfile();
00687
00694 virtual SegmentProfile& operator=(const SegmentProfile& clone);
00695
00696 protected:
00697
00698 std::vector<std::vector<MSK_ELM_TYPE> > elements_;
00700 const ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>* elmProfile_;
00704 std::map<MSK_ELM_TYPE, MSK_ELM_TYPE> mask_;
00706 std::vector<unsigned int> idxMask_;
00708 std::map<MSK_ELM_TYPE, unsigned int> age_;
00710 std::map<MSK_ELM_TYPE, bool> inactive_;
00712 std::map<MSK_ELM_TYPE, MSK_ELM_TYPE> simSeg_;
00715 unsigned int nSeg_;
00716 unsigned int nChl_;
00717 MSK_ELM_TYPE nextSegId_;
00718
00719 Matrix neighbors_;
00722
00723 double** intensity_;
00724 double* locationX_;
00725 double* locationY_;
00726 unsigned int* nElements_;
00727 unsigned int* nPixels_;
00729 private:
00730
00740 void attachIds(const SegmentProfile<IMG_ELM_TYPE,MSK_ELM_TYPE>& last);
00741
00750 void attach(const SegmentProfile<IMG_ELM_TYPE,MSK_ELM_TYPE>& last);
00751
00755 void cleanup();
00756 };
00757
00758
00759
00760
00761 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
00762 unsigned int Profile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
00763 idx(const MSK_ELM_TYPE& id) const
00764 {
00765 typename std::map<MSK_ELM_TYPE, unsigned int>::const_iterator elmIdx_p;
00766
00767 elmIdx_p = indices_.find(id);
00768 if (elmIdx_p == indices_.end())
00769 throw ncutException("element not found");
00770
00771 return elmIdx_p->second;
00772 }
00773
00774 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
00775 Profile<IMG_ELM_TYPE, MSK_ELM_TYPE>::Profile(const Setting* setting)
00776 : similarity_(), indices_(), identifiers_(), setting_(setting)
00777 {
00778 return;
00779 }
00780
00781 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
00782 Profile<IMG_ELM_TYPE, MSK_ELM_TYPE>::Profile(const Profile& clone)
00783 : similarity_(), indices_(), identifiers_(), setting_(clone.setting_)
00784 {
00785 *this = clone;
00786 }
00787
00788 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
00789 Profile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
00790 Profile(const Matrix& similarity,
00791 const std::map<MSK_ELM_TYPE, unsigned int>& indices,
00792 const std::vector<MSK_ELM_TYPE>& identifiers,
00793 const Setting* setting)
00794 : similarity_(similarity), indices_(indices), identifiers_(identifiers),
00795 setting_(setting)
00796 {
00797 if (similarity_.dim() == 0 ||
00798 similarity_.dim() != indices_.size() ||
00799 similarity_.dim() != identifiers_.size()) {
00800
00801 throw ncutException("dimension = 0 or dimension mismatch");
00802 }
00803 }
00804
00805 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
00806 Profile<IMG_ELM_TYPE, MSK_ELM_TYPE>::~Profile()
00807 {
00808 return;
00809 }
00810
00811 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
00812 Profile<IMG_ELM_TYPE,MSK_ELM_TYPE>& Profile<IMG_ELM_TYPE,MSK_ELM_TYPE>::
00813 operator=(const Profile& clone)
00814 {
00815
00816 if (this == &clone) return *this;
00817
00818
00819 if (!indices_.empty()) indices_.clear();
00820 if (!identifiers_.empty()) identifiers_.clear();
00821
00822
00823 this->similarity_.changeDimension(clone.similarity_.dim());
00824 this->similarity_ = clone.similarity_;
00825
00826
00827 this->indices_ = clone.indices_;
00828 this->identifiers_ = clone.identifiers_;
00829
00830
00831 this->setting_ = clone.setting_;
00832
00833 return *this;
00834 }
00835
00836 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
00837 unsigned int ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
00838 frameOf(unsigned int idx) const
00839 {
00840 for (unsigned int frm=0; frm<maxElmIdx_.size(); frm++)
00841 if (maxElmIdx_[frm] >= idx) return frm;
00842
00843 throw ncutException("element index out of bounds");
00844 }
00845
00846 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
00847 void ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::readaptBase()
00848 {
00849 #ifdef debugout
00850 std::cout << "readaptBase" << std::endl << std::flush;
00851 #endif
00852 nElm_front_ = maxElmIdx_.front() + 1;
00853 nElm_old_ = nElm_;
00854
00855
00856 typename std::map<MSK_ELM_TYPE, unsigned int>::iterator indicesItr;
00857 typename std::map<MSK_ELM_TYPE, unsigned int>::iterator eraseItr;
00858 for (indicesItr = this->indices_.begin();
00859 indicesItr != this->indices_.end(); ) {
00860 if (indicesItr->second < nElm_front_) {
00861 eraseItr = indicesItr++;
00862 this->indices_.erase(eraseItr);
00863 }
00864 else {
00865 indicesItr->second -= nElm_front_;
00866 indicesItr++;
00867 }
00868 }
00869
00870
00871 this->identifiers_.resize(this->indices_.size());
00872 for (indicesItr = this->indices_.begin();
00873 indicesItr != this->indices_.end();
00874 indicesItr++) {
00875 this->identifiers_[indicesItr->second] = indicesItr->first;
00876 }
00877
00878
00879 maxElmIdx_.erase(maxElmIdx_.begin());
00880 typename std::vector<unsigned int>::iterator maxElmIdxItr;
00881 for (maxElmIdxItr=maxElmIdx_.begin();
00882 maxElmIdxItr!=maxElmIdx_.end();
00883 maxElmIdxItr++) {
00884 (*maxElmIdxItr) -= nElm_front_;
00885 }
00886
00887
00888 typename std::vector<std::vector<unsigned int> >::iterator pxlsItr =
00889 pxls_.begin();
00890 for (unsigned int elm=0; elm <nElm_front_; elm++)
00891 pxlsItr++;
00892 pxls_.erase(pxls_.begin(), pxlsItr);
00893
00894 }
00895
00896 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
00897 void ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
00898 chartElements(unsigned int frmNo)
00899 {
00900 #ifdef debugout
00901 std::cout << "chartElements" << std::endl << std::flush;
00902 #endif
00903
00904 if (frmNo >= maskSeq_.length())
00905 throw ncutException("invalid frame number, corrupt mask sequence?");
00906
00907 Image<MSK_ELM_TYPE>* maskImg = maskSeq_[frmNo];
00908 unsigned int maskPxlCount = width_*height_;
00909 MSK_ELM_TYPE currentElementName;
00910 unsigned int startNumber;
00911 if (maxElmIdx_.empty())
00912 startNumber = 0;
00913 else
00914 startNumber = maxElmIdx_.back()+1;
00915 unsigned int elmNumber = startNumber;
00916
00917 typename std::map<MSK_ELM_TYPE, unsigned int>::iterator firstOccurrence;
00918
00919 for (unsigned int pxl=0; pxl<maskPxlCount; pxl++) {
00920 currentElementName = maskImg->val(pxl);
00921 firstOccurrence = this->indices_.find(currentElementName);
00922 if (firstOccurrence == this->indices_.end()) {
00923
00924 this->identifiers_.push_back(currentElementName);
00925 this->indices_[currentElementName] = elmNumber++;
00926 this->pxls_.push_back(std::vector<unsigned int>());
00927 this->pxls_.back().push_back(pxl);
00928 }
00929 else if (firstOccurrence->second < startNumber) {
00930
00931 throw ncutException("duplicate element name");
00932 }
00933 else {
00934 this->pxls_[firstOccurrence->second].push_back(pxl);
00935 }
00936 }
00937 maxElmIdx_.push_back(elmNumber-1);
00938
00939 nElm_ = elmNumber;
00940
00941 if (nElm_ <= 0)
00942 throw ncutException("no elements found");
00943 }
00944
00945 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
00946 void ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::initInfoArrays()
00947 {
00948 #ifdef debugout
00949 std::cout << "initInfoArrays" << std::endl << std::flush;
00950 #endif
00951
00952
00953 hue_ = new double [nElm_];
00954 hueVariation_ = new double [nElm_];
00955 hueVecX_ = new double [nElm_];
00956 hueVecY_ = new double [nElm_];
00957 lightness_ = new double [nElm_];
00958 litVariation_ = new double [nElm_];
00959 saturation_ = new double [nElm_];
00960 satVariation_ = new double [nElm_];
00961 intensity_ = new double* [nElm_];
00962 for (unsigned int elm=0; elm<nElm_; elm++)
00963 intensity_[elm] = new double[frameSeq_.nChl()];
00964 locationX_ = new double [nElm_];
00965 locationY_ = new double [nElm_];
00966 locationZ_ = new double [nElm_];
00967 nPixels_ = new unsigned int [nElm_];
00968
00969
00970 for (unsigned int elm=0; elm<nElm_; elm++) {
00971 hue_[elm] = 0.0;
00972 hueVariation_[elm] = 0.0;
00973 hueVecX_[elm] = 0.0;
00974 hueVecY_[elm] = 0.0;
00975 lightness_[elm] = 0.0;
00976 litVariation_[elm] = 0.0;
00977 saturation_[elm] = 0.0;
00978 satVariation_[elm] = 0.0;
00979 for (unsigned int chl=0; chl<frameSeq_.nChl(); chl++)
00980 intensity_[elm][chl] = 0.0;
00981 locationX_[elm] = 0.0;
00982 locationY_[elm] = 0.0;
00983 locationZ_[elm] = 0.0;
00984 nPixels_[elm] = 0;
00985 }
00986 }
00987
00988 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
00989 void ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::initInfo()
00990 {
00991 #ifdef debugout
00992 std::cout << "initInfo" << std::endl << std::flush;
00993 #endif
00994
00995
00996 this->similarity_.changeDimension(nElm_);
00997 this->neighbors_.changeDimension(nElm_);
00998
00999
01000 initInfoArrays();
01001 }
01002
01003 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
01004 void ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::readaptInfo()
01005 {
01006 #ifdef debugout
01007 std::cout << "readaptInfo" << std::endl << std::flush;
01008 #endif
01009
01010
01011 Matrix firstFrameSim, otherFramesSim;
01012 Matrix firstFrameNeighb, otherFramesNeighb;
01013 Matrix firstFrameBordrs, otherFramesBordrs;
01014
01015 std::vector<unsigned int> indexMap;
01016 std::list<unsigned int> firstFrameIndices;
01017
01018 for (unsigned int elm=0; elm<nElm_front_; elm++)
01019 firstFrameIndices.push_back(elm);
01020
01021 this->similarity_.split(firstFrameSim, otherFramesSim, indexMap,
01022 firstFrameIndices);
01023 this->neighbors_.split(firstFrameNeighb, otherFramesNeighb, indexMap,
01024 firstFrameIndices);
01025
01026 this->similarity_.changeDimension(nElm_);
01027 this->neighbors_.changeDimension(nElm_);
01028
01029 this->similarity_.inject(otherFramesSim, 0);
01030 this->neighbors_.inject(otherFramesNeighb, 0);
01031
01032
01033 double* hue_old = hue_;
01034 double* hueVariation_old = hueVariation_;
01035 double* hueVecX_old = hueVecX_;
01036 double* hueVecY_old = hueVecY_;
01037 double* lightness_old = lightness_;
01038 double* litVariation_old = litVariation_;
01039 double* saturation_old = saturation_;
01040 double* satVariation_old = satVariation_;
01041 double** intensity_old = intensity_;
01042 double* locationX_old = locationX_;
01043 double* locationY_old = locationY_;
01044 double* locationZ_old = locationZ_;
01045 unsigned int* nPixels_old = nPixels_;
01046
01047
01048 initInfoArrays();
01049
01050
01051
01052 unsigned int old_i, new_i;
01053 for (old_i=nElm_front_, new_i=0; old_i<nElm_old_; old_i++, new_i++) {
01054 hue_[new_i] = hue_old[old_i];
01055 hueVariation_[new_i] = hueVariation_old[old_i];
01056 hueVecX_[new_i] = hueVecX_old[old_i];
01057 hueVecY_[new_i] = hueVecY_old[old_i];
01058 lightness_[new_i] = lightness_old[old_i];
01059 litVariation_[new_i] = litVariation_old[old_i];
01060 saturation_[new_i] = saturation_old[old_i];
01061 satVariation_[new_i] = satVariation_old[old_i];
01062 for (unsigned int chl=0; chl<frameSeq_.nChl(); chl++) {
01063 intensity_[new_i][chl] = intensity_old[old_i][chl];
01064 }
01065 locationX_[new_i] = locationX_old[old_i];
01066 locationY_[new_i] = locationY_old[old_i];
01067 locationZ_[new_i] = locationZ_old[old_i] - 1.0;
01068 nPixels_[new_i] = nPixels_old[old_i];
01069 }
01070
01071
01072 if (hue_old != NULL) {
01073 delete [] hue_old;
01074 hue_old = NULL;
01075 }
01076 if (hueVariation_old != NULL) {
01077 delete [] hueVariation_old;
01078 hueVariation_old = NULL;
01079 }
01080 if (hueVecX_old != NULL) {
01081 delete [] hueVecX_old;
01082 hueVecX_old = NULL;
01083 }
01084 if (hueVecY_old != NULL) {
01085 delete [] hueVecY_old;
01086 hueVecY_old = NULL;
01087 }
01088 if (lightness_old != NULL) {
01089 delete [] lightness_old;
01090 lightness_old = NULL;
01091 }
01092 if (litVariation_old != NULL) {
01093 delete [] litVariation_old;
01094 litVariation_old = NULL;
01095 }
01096 if (saturation_old != NULL) {
01097 delete [] saturation_old;
01098 saturation_old = NULL;
01099 }
01100 if (satVariation_old != NULL) {
01101 delete [] satVariation_old;
01102 satVariation_old = NULL;
01103 }
01104 if (intensity_old != NULL) {
01105 for (unsigned int elm=0; elm < nElm_old_; elm++) {
01106 if (intensity_old[elm] != NULL) delete [] intensity_old[elm];
01107 }
01108 delete [] intensity_old;
01109 intensity_old = NULL;
01110 }
01111 if (locationX_old != NULL) {
01112 delete [] locationX_old;
01113 locationX_old = NULL;
01114 }
01115 if (locationY_old != NULL) {
01116 delete [] locationY_old;
01117 locationY_old = NULL;
01118 }
01119 if (locationZ_old != NULL) {
01120 delete [] locationZ_old;
01121 locationZ_old = NULL;
01122 }
01123 if (nPixels_old != NULL) {
01124 delete [] nPixels_old;
01125 nPixels_old = NULL;
01126 }
01127 }
01128
01129 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
01130 void ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
01131 getInfo(unsigned int from, unsigned int to)
01132 {
01133 #ifdef debugout
01134 std::cout << "getInfo" << std::endl << std::flush;
01135 #endif
01136
01137 if (to >= frameSeq_.length() || from > to)
01138 throw ncutException("invalid frame range, corrupt sequences?");
01139
01140 unsigned int varBins = this->setting_->variation_bins;
01141 double satMult = this->setting_->saturation_multiplier;
01142
01143 if (varBins < 3 || varBins > IMAGE_DATATYPE_MAXVAL)
01144 throw ncutException(
01145 "invalid number of bins, must be in [3,IMAGE_DATATYPE_MAXVAL]" );
01146 if (satMult < 0.0)
01147 throw ncutException("invalid saturation multiplier, must be >= 0");
01148
01149
01150 const double binFac = normFac * varBins;
01151
01152
01153 const unsigned int pxlCount = width_*height_;
01154
01155
01156 const unsigned int cmpCount = pxlCount*frameSeq_.nChl();
01157
01158
01159 std::deque<pxlInfo_t> pxlHist;
01160 unsigned int pxlHistSize = (3*width_)+3;
01161 unsigned int pxlHistOffset = width_+1;
01162 if (pxlHistOffset >= pxlHistSize)
01163 throw ncutException(
01164 "invalid pixel history offset specified (offset >= size)" );
01165
01166
01167 double** satBuffer = new double* [nElm_];
01168 for (unsigned int i=0; i<nElm_; i++) {
01169 satBuffer[i] = new double[varBins];
01170 for (unsigned int j=0; j<varBins; j++) {
01171 satBuffer[i][j] = 0.0;
01172 }
01173 }
01174
01175
01176 double** litBuffer = new double* [nElm_];
01177 for (unsigned int i=0; i<nElm_; i++) {
01178 litBuffer[i] = new double[varBins];
01179 for (unsigned int j=0; j<varBins; j++) {
01180 litBuffer[i][j] = 0.0;
01181 }
01182 }
01183
01184
01185 double intVal [frameSeq_.nChl()];
01186
01187
01188 unsigned int elm;
01189 unsigned int chl;
01190 unsigned int cmp;
01191 unsigned int pxl;
01192 unsigned int col;
01193 unsigned int row;
01194 unsigned int histPxl;
01195 unsigned int histCol;
01196 unsigned int histRow;
01197
01198 for (unsigned int frm=from; frm<=to; frm++) {
01199
01200
01201 chl = 0;
01202 cmp = 0;
01203 pxl = 0;
01204 col = 0;
01205 row = 0;
01206 histPxl = 0;
01207 histCol = 0;
01208 histRow = 0;
01209
01210
01211 pxlHist.clear();
01212
01213
01214 Image<MSK_ELM_TYPE>* mask = maskSeq_[frm];
01215 Image<IMG_ELM_TYPE>* image = frameSeq_[frm];
01216
01217 while (histPxl<pxlCount) {
01218
01219
01220 if (cmp < cmpCount) {
01221
01222
01223 if (chl == 0) {
01224 typename std::map<MSK_ELM_TYPE, unsigned int>::iterator currentElmItr
01225 = this->indices_.find(mask->val(pxl));
01226 if (currentElmItr == this->indices_.end())
01227 throw ncutCritical(
01228 "data corrupt: unknown element, index map is incomplete" );
01229 elm = currentElmItr->second;
01230
01231 if (elm >= nElm_)
01232 throw ncutCritical(
01233 "data corrupt: invalid element index in index map" );
01234 }
01235
01236
01237 IMG_ELM_TYPE curInt = image->val(cmp);
01238 intVal[chl] = curInt;
01239 intensity_[elm][chl] += curInt;
01240
01241
01242 if (chl == frameSeq_.nChl()-1) {
01243
01244
01245 locationX_[elm] += col;
01246 locationY_[elm] += row;
01247 locationZ_[elm] += frm;
01248
01249
01250 nPixels_[elm]++;
01251
01252
01253 double pxlLightness = (0.2126*intVal[0]) +
01254 (0.7152*intVal[1]) +
01255 (0.0722*intVal[2]);
01256 lightness_[elm] += pxlLightness;
01257
01258
01259 double maxIntVal, minIntVal;
01260 maxIntVal = minIntVal = intVal[0];
01261 for (unsigned int i=1; i<frameSeq_.nChl(); i++) {
01262 if (intVal[i] > maxIntVal) maxIntVal = intVal[i];
01263 if (intVal[i] < minIntVal) minIntVal = intVal[i];
01264 }
01265 double pxlSaturation = maxIntVal - minIntVal;
01266 pxlSaturation *= satMult;
01267 if (pxlSaturation > 255.0) pxlSaturation = 255.0;
01268
01269 saturation_[elm] += pxlSaturation;
01270
01271
01272 double pxlHue, pxlHueVecX, pxlHueVecY;
01273 if (pxlSaturation != 0.0) {
01274 pxlHue = acos((intVal[0] - (intVal[1]/2.0) - (intVal[2]/2.0)) /
01275 sqrt((intVal[0]*intVal[0]) + (intVal[1]*intVal[1]) +
01276 (intVal[2]*intVal[2]) - (intVal[0]*intVal[1]) -
01277 (intVal[0]*intVal[2]) - (intVal[2]*intVal[1])));
01278 pxlHue = (intVal[2] > intVal[1]) ? ((2*PI) - pxlHue) : pxlHue;
01279
01280 pxlHueVecX = pxlSaturation * cos(pxlHue);
01281 pxlHueVecY = pxlSaturation * sin(pxlHue);
01282 }
01283 else {
01284 pxlHueVecX = 0.0;
01285 pxlHueVecY = 0.0;
01286 }
01287 hueVecX_[elm] += pxlHueVecX;
01288 hueVecY_[elm] += pxlHueVecY;
01289
01290
01291
01292 double binPos = (pxlSaturation * binFac) - 0.5;
01293 int binNum;
01294 if (binPos < 0.0) satBuffer[elm][0]++;
01295 else if (binPos >= (double) (varBins-1))
01296 satBuffer[elm][varBins-1]++;
01297 else {
01298 binNum = (int) binPos;
01299 binPos -= binNum;
01300 satBuffer[elm][binNum+1]+=binPos;
01301 satBuffer[elm][binNum]+=1.0-binPos;
01302 }
01303 binPos = (pxlLightness * binFac) - 0.5;
01304 if (binPos < 0.0) litBuffer[elm][0]++;
01305 else if (binPos >= (double) (varBins-1))
01306 litBuffer[elm][varBins-1]++;
01307 else {
01308 binNum = (int) binPos;
01309 binPos -= binNum;
01310 litBuffer[elm][binNum+1]+=binPos;
01311 litBuffer[elm][binNum]+=1.0-binPos;
01312 }
01313
01314
01315 pxlInfo_t pxlInfo;
01316 pxlInfo.hueVecX = pxlHueVecX;
01317 pxlInfo.hueVecY = pxlHueVecY;
01318 pxlInfo.lit = pxlLightness;
01319 pxlInfo.elm = elm;
01320
01321 if (pxl >= pxlHistSize) pxlHist.pop_back();
01322 pxlHist.push_front(pxlInfo);
01323 }
01324 }
01325 else {
01326
01327
01328
01329 if (chl == frameSeq_.nChl()-1) {
01330 pxlInfo_t pxlInfo;
01331 pxlInfo.hueVecX = 0.0;
01332 pxlInfo.hueVecY = 0.0;
01333 pxlInfo.lit = 0.0;
01334 pxlInfo.elm = 0;
01335
01336 if (pxl >= pxlHistSize) pxlHist.pop_back();
01337 pxlHist.push_front(pxlInfo);
01338 }
01339 }
01340
01341
01342
01343 if (chl == frameSeq_.nChl()-1) {
01344
01345
01346 if (histPxl >= width_) {
01347 unsigned int thisElm = pxlHist[pxlHistOffset].elm;
01348 unsigned int neighborElm = pxlHist[pxlHistOffset+width_].elm;
01349 if (neighborElm != thisElm) {
01350 neighbors_.val(neighborElm,thisElm) =
01351 ++neighbors_.val(thisElm,neighborElm);
01352 }
01353 }
01354 if (histCol >= 1) {
01355 unsigned int thisElm = pxlHist[pxlHistOffset].elm;
01356 unsigned int neighborElm = pxlHist[pxlHistOffset+1].elm;
01357 if (neighborElm != thisElm) {
01358 neighbors_.val(neighborElm,thisElm) =
01359 ++neighbors_.val(thisElm,neighborElm);
01360 }
01361 }
01362 if (frm >= 1 && pxl>=pxlHistOffset) {
01363 unsigned int thisElm = pxlHist[pxlHistOffset].elm;
01364 typename std::map<MSK_ELM_TYPE, unsigned int>::iterator lastFrmElmItr
01365 = this->indices_.find(maskSeq_[frm-1]->val(histPxl));
01366 if (lastFrmElmItr == this->indices_.end())
01367 throw ncutException(
01368 "data corrupt: invalid element index in index map" );
01369 unsigned int lastFrmElm = lastFrmElmItr->second;
01370 if (lastFrmElm != thisElm) {
01371 if (neighbors_.val(lastFrmElm,thisElm) < 0.5)
01372 neighbors_.val(lastFrmElm,thisElm) =
01373 neighbors_.val(thisElm,lastFrmElm) = 1.0;
01374 }
01375 }
01376 }
01377
01378
01379 if (chl == frameSeq_.nChl()-1) {
01380 if (pxl>=pxlHistOffset) {
01381 histPxl++;
01382 if (++histCol >= width_) {
01383 histCol = 0;
01384 histRow++;
01385 }
01386 }
01387
01388 pxl++;
01389
01390 if (++col >= width_) {
01391 col = 0;
01392 row++;
01393 }
01394 }
01395
01396 if (++chl >= frameSeq_.nChl())
01397 chl = 0;
01398
01399 cmp++;
01400 }
01401 }
01402
01403
01404 for (unsigned int elm=minIdx(from); elm<=maxIdx(to); elm++) {
01405
01406
01407 if (nPixels_[elm]==0) {
01408
01409 hue_[elm] = 0.0;
01410 hueVecX_[elm] = 0.0;
01411 hueVecY_[elm] = 0.0;
01412 hueVariation_[elm] = 0.0;
01413 lightness_[elm] = 0.0;
01414 saturation_[elm] = 0.0;
01415 litVariation_[elm] = 0.0;
01416 satVariation_[elm] = 0.0;
01417 for (unsigned int chl=0; chl<frameSeq_.nChl(); chl++)
01418 intensity_[elm][chl] = 0;
01419 locationX_[elm] = 0;
01420 locationY_[elm] = 0;
01421 locationZ_[elm] = 0;
01422
01423 continue;
01424 }
01425
01426
01427 if (hueVecX_[elm] > 0.0) {
01428 if (hueVecY_[elm] >= 0.0)
01429 hue_[elm] = atan(hueVecY_[elm] / hueVecX_[elm]);
01430 else
01431 hue_[elm] = atan(hueVecY_[elm] / hueVecX_[elm]) + (2*PI);
01432 }
01433 else if (hueVecX_[elm] < 0.0) {
01434 hue_[elm] = atan(hueVecY_[elm] / hueVecX_[elm]) + PI;
01435 }
01436 else {
01437 if (hueVecY_[elm] > 0.0)
01438 hue_[elm] = PI*0.5;
01439 else if (hueVecY_[elm] < 0.0)
01440 hue_[elm] = PI*1.5;
01441 else
01442 hue_[elm] = -1.0;
01443 }
01444 if (hue_[elm] > 0.0) hue_[elm] *= 180.0 / PI;
01445
01446
01447 if (saturation_[elm] < 0.01)
01448
01449 hueVariation_[elm] = 1.0;
01450 else
01451
01452 hueVariation_[elm] = sqrt((hueVecX_[elm]*hueVecX_[elm]) +
01453 (hueVecY_[elm]*hueVecY_[elm]) )
01454 / saturation_[elm];
01455 hueVariation_[elm] = 1.0 - hueVariation_[elm];
01456 if (hueVariation_[elm] < 0.0) hueVariation_[elm] = 0.0;
01457
01458
01459 hueVecX_[elm] /= nPixels_[elm];
01460 hueVecY_[elm] /= nPixels_[elm];
01461
01462
01463 lightness_[elm] /= nPixels_[elm];
01464 saturation_[elm] /= nPixels_[elm];
01465
01466
01467 for (unsigned int bin=0; bin<varBins; bin++) {
01468 if (litBuffer[elm][bin] <= 0.0) continue;
01469 litVariation_[elm] += fabs( (((double)bin)/binFac) - lightness_[elm])
01470 * litBuffer[elm][bin];
01471 }
01472 litVariation_[elm] /= nPixels_[elm];
01473 litVariation_[elm] /= IMAGE_DATATYPE_MAXVAL;
01474
01475 for (unsigned int bin=0; bin<varBins; bin++) {
01476 if (satBuffer[elm][bin] <= 0.0) continue;
01477 satVariation_[elm] += fabs( (((double)bin)/binFac) - saturation_[elm])
01478 * satBuffer[elm][bin];
01479 }
01480 satVariation_[elm] /= nPixels_[elm];
01481 satVariation_[elm] /= IMAGE_DATATYPE_MAXVAL;
01482
01483
01484 for (unsigned int chl=0; chl<frameSeq_.nChl(); chl++) {
01485 intensity_[elm][chl] = intensity_[elm][chl]/nPixels_[elm];
01486 }
01487
01488
01489 locationX_[elm] = locationX_[elm]/nPixels_[elm];
01490 locationY_[elm] = locationY_[elm]/nPixels_[elm];
01491 locationZ_[elm] = locationZ_[elm]/nPixels_[elm];
01492 }
01493
01494
01495 if (litBuffer != NULL) {
01496 for (unsigned int elm=0; elm<nElm_; elm++) {
01497 if (litBuffer[elm] != NULL)
01498 delete [] litBuffer[elm];
01499 }
01500 delete [] litBuffer;
01501 litBuffer = NULL;
01502 }
01503 if (satBuffer != NULL) {
01504 for (unsigned int elm=0; elm<nElm_; elm++) {
01505 if (satBuffer[elm] != NULL)
01506 delete [] satBuffer[elm];
01507 }
01508 delete [] satBuffer;
01509 satBuffer = NULL;
01510 }
01511 }
01512
01513 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
01514 void ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
01515 calcSim(unsigned int from, unsigned int to)
01516 {
01517 #ifdef debugout
01518 std::cout << "calcSim" << std::endl << std::flush;
01519 #endif
01520
01521 if (to >= frameSeq_.length() || from > to)
01522 throw ncutException("invalid frame range, corrupt sequences?");
01523
01524
01525 double intWeight = this->setting_->intensity_weight;
01526
01527 double posWeight = this->setting_->position_weight;
01528
01529 double falloff_int = this->setting_->intensity_falloff;
01530
01531 double falloff_pos = this->setting_->distance_falloff;
01532
01533 bool stretchSat = (this->setting_->stretch_saturation == 1);
01534
01535 if (intWeight < 0.0 || intWeight > 1.0)
01536 throw ncutException("invalid intensity weight, must be in [0,1]");
01537 if (posWeight < 0.0 || posWeight > 1.0)
01538 throw ncutException("invalid position weight, must be in [0,1]");
01539 if (falloff_int < 0.000001)
01540 throw ncutException("invalid intensity falloff, must be >= 0.000001");
01541 if (falloff_pos < 0.000001)
01542 throw ncutException("invalid distance falloff, must be >= 0.000001");
01543
01544
01545 double point0_x = this->setting_->stretch_point0_x;
01546 double point0_y = this->setting_->stretch_point0_y;
01547 double point1_x = this->setting_->stretch_point1_x;
01548 double point1_y = this->setting_->stretch_point1_y;
01549 double point2_x = this->setting_->stretch_point2_x;
01550 double point2_y = this->setting_->stretch_point2_y;
01551 double point3_x = this->setting_->stretch_point3_x;
01552 double point3_y = this->setting_->stretch_point3_y;
01553
01554 if (point0_x < 0.0 || point0_x > IMAGE_DATATYPE_MAXVAL)
01555 throw ncutException("invalid stretch point 0 coordinate x");
01556 if (point0_y < 0.0 || point0_y > IMAGE_DATATYPE_MAXVAL)
01557 throw ncutException("invalid stretch point 0 coordinate y");
01558 if (point1_x <= point0_x || point1_x > IMAGE_DATATYPE_MAXVAL)
01559 throw ncutException("invalid stretch point 1 coordinate x");
01560 if (point1_y < 0.0 || point1_y > IMAGE_DATATYPE_MAXVAL)
01561 throw ncutException("invalid stretch point 1 coordinate y");
01562 if (point2_x <= point1_x || point2_x > IMAGE_DATATYPE_MAXVAL)
01563 throw ncutException("invalid stretch point 2 coordinate x");
01564 if (point2_y < 0.0 || point2_y > IMAGE_DATATYPE_MAXVAL)
01565 throw ncutException("invalid stretch point 2 coordinate y");
01566 if (point3_x <= point2_x || point3_x > IMAGE_DATATYPE_MAXVAL)
01567 throw ncutException("invalid stretch point 3 coordinate x");
01568 if (point3_y < 0.0 || point3_y > IMAGE_DATATYPE_MAXVAL)
01569 throw ncutException("invalid stretch point 3 coordinate y");
01570
01571
01572 const double perctFac = 100.0*normFac;
01573
01574
01575 double intSim;
01576 double posSim;
01577
01578 unsigned int firstIdx = minIdx(from);
01579 unsigned int lastIdx = maxIdx(to);
01580
01581 for (unsigned int i=0; i<=lastIdx; i++) {
01582 this->similarity_.val(i,i) = 0.0;
01583 for (unsigned int j=(i<firstIdx)?firstIdx:(i+1); j<=lastIdx; j++) {
01584
01585 if (intWeight > 0.0) {
01586 double hueVecX_i, hueVecY_i, hueVecX_j, hueVecY_j;
01587 if (stretchSat) {
01588 double sat_i = saturation_[i];
01589 double sat_j = saturation_[j];
01590
01591 double interPos;
01592 double stretch_i;
01593 if (sat_i <= point1_x) {
01594 interPos = (sat_i-point0_x) / (point1_x-point0_x);
01595 stretch_i = (interPos*point1_y) + (1-interPos)*point0_y;
01596 }
01597 else if (sat_i <= point2_x) {
01598 interPos = (sat_i-point1_x) / (point2_x-point1_x);
01599 stretch_i = (interPos*point2_y) + (1-interPos)*point1_y;
01600 }
01601 else {
01602 interPos = (sat_i-point2_x) / (point3_x-point2_x);
01603 stretch_i = (interPos*point3_y) + (1-interPos)*point2_y;
01604 }
01605
01606 double stretch_j;
01607 if (sat_j <= point1_x) {
01608 interPos = (sat_j-point0_x) / (point1_x-point0_x);
01609 stretch_j = (interPos*point1_y) + (1-interPos)*point0_y;
01610 }
01611 else if (sat_j <= point2_x) {
01612 interPos = (sat_j-point1_x) / (point2_x-point1_x);
01613 stretch_j = (interPos*point2_y) + (1-interPos)*point1_y;
01614 }
01615 else {
01616 interPos = (sat_j-point2_x) / (point3_x-point2_x);
01617 stretch_j = (interPos*point3_y) + (1-interPos)*point2_y;
01618 }
01619
01620 if (sat_i > 0.0) {
01621 hueVecX_i = hueVecX_[i] * (stretch_i/sat_i);
01622 hueVecY_i = hueVecY_[i] * (stretch_i/sat_i);
01623 }
01624 else {
01625 hueVecX_i = 0.0;
01626 hueVecY_i = 0.0;
01627 }
01628 if (sat_j > 0.0) {
01629 hueVecX_j = hueVecX_[j] * (stretch_j/sat_j);
01630 hueVecY_j = hueVecY_[j] * (stretch_j/sat_j);
01631 }
01632 else {
01633 hueVecX_j = 0.0;
01634 hueVecY_j = 0.0;
01635 }
01636 }
01637 else {
01638 hueVecX_i = hueVecX_[i];
01639 hueVecY_i = hueVecY_[i];
01640 hueVecX_j = hueVecX_[j];
01641 hueVecY_j = hueVecY_[j];
01642 }
01643
01644 double hueSatDifX = fabs(hueVecX_i - hueVecX_j);
01645 hueSatDifX *= perctFac;
01646
01647 double hueSatDifY = fabs(hueVecY_i - hueVecY_j);
01648 hueSatDifY *= perctFac;
01649
01650 double litDif = fabs(lightness_[i] - lightness_[j]);
01651 litDif *= perctFac;
01652
01653 double distance = sqrt(hueSatDifX*hueSatDifX +
01654 hueSatDifY*hueSatDifY +
01655 litDif*litDif);
01656 intSim = sqrt(200.0*200.0 +
01657 200.0*200.0 +
01658 100.0*100.0) - distance;
01659
01660
01661
01662 if (intSim < 0.0) intSim = 0.0;
01663 intSim = (1.0-intWeight)+(intSim*intWeight);
01664 }
01665 else
01666 intSim = 1.0;
01667
01668
01669 if (posWeight > 0.0) {
01670 if (neighbors_.val(i,j) > 0.5)
01671 posSim = 1.0;
01672 else {
01673 double sigma = falloff_pos * ((width_ < height_) ? width_:height_);
01674 double distance = sqrt((fabs(locationX_[i] - locationX_[j])*
01675 fabs(locationX_[i] - locationX_[j])) +
01676 (fabs(locationY_[i] - locationY_[j])*
01677 fabs(locationY_[i] - locationY_[j])) +
01678 (fabs(locationZ_[i] - locationZ_[j])*
01679 fabs(locationZ_[i] - locationZ_[j])));
01680 posSim = exp(-1.0 * pow(distance/sigma, 2.0));
01681 }
01682 if (posSim < 0.0) posSim = 0.0;
01683 posSim = (1.0-posWeight)+(posSim*posWeight);
01684 }
01685 else
01686 posSim = 1.0;
01687
01688 this->similarity_.val(i,j) = this->similarity_.val(j,i) = intSim * posSim;
01689 }
01690 }
01691 }
01692
01693 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
01694 int ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
01695 push(const Image<IMG_ELM_TYPE>* frame, const Image<MSK_ELM_TYPE>* mask)
01696 {
01697 if (frame == NULL || mask == NULL ||
01698 frame->height() != mask->height() ||
01699 frame->width() != mask->width() ||
01700 mask->nChl() != 1 ||
01701 frame->nChl() != 3)
01702 throw ncutException("wrong image format");
01703
01704 if (frameSeq_.length() == 0) {
01705 this->width_ = frame->width();
01706 this->height_ = frame->height();
01707 this->nChl_ = frame->nChl();
01708 }
01709
01710 try {
01711 if (frameSeq_.length() < wSize_) {
01712 frameSeq_.pushFrame(*frame);
01713 maskSeq_.pushFrame(*mask);
01714 }
01715 else {
01716 frameSeq_.pushThrough(*frame);
01717 maskSeq_.pushThrough(*mask);
01718 }
01719 }
01720 catch(ncutException& e) {
01721 if (frameSeq_.length() > maskSeq_.length()) frameSeq_.popFrame();
01722 throw ncutException(
01723 (std::string("wrong image format : ")+e.what()).c_str() );
01724 }
01725
01726 if (frameSeq_.length() < wSize_) {
01727 chartElements(maskSeq_.length()-1);
01728
01729 return 1;
01730 }
01731 else {
01732 if (this->similarity_.dim() == 0) {
01733 chartElements(maskSeq_.length()-1);
01734 initInfo();
01735 getInfo(0, frameSeq_.length()-1);
01736 calcSim(0, frameSeq_.length()-1);
01737 }
01738 else {
01739 readaptBase();
01740 chartElements(maskSeq_.length()-1);
01741 readaptInfo();
01742 getInfo(frameSeq_.length()-1, frameSeq_.length()-1);
01743 calcSim(frameSeq_.length()-1, frameSeq_.length()-1);
01744 }
01745
01746 return 0;
01747 }
01748 }
01749
01750 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
01751 ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
01752 ElementProfile(unsigned int wSize, const Setting* setting)
01753 : Profile<IMG_ELM_TYPE, MSK_ELM_TYPE>(setting),
01754 normFac(1.0/IMAGE_DATATYPE_MAXVAL),
01755 frameSeq_(), maskSeq_(),
01756 wSize_(wSize), width_(0), height_(0), nChl_(0),
01757 nElm_(0), pxls_(), maxElmIdx_(), neighbors_(),
01758 hue_(NULL), hueVariation_(NULL), hueVecX_(NULL), hueVecY_(NULL),
01759 lightness_(NULL), litVariation_(NULL), saturation_(NULL),
01760 satVariation_(NULL), intensity_(NULL), locationX_(NULL), locationY_(NULL),
01761 locationZ_(NULL), nPixels_(NULL), nElm_old_(0), nElm_front_(0)
01762
01763 {
01764 if (wSize_ < 2 || wSize_ % 2 != 1)
01765 throw ncutException("invalid window size, wSize must be >= 2 and odd");
01766
01767 return;
01768 }
01769
01770 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
01771 ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
01772 ElementProfile(const ElementProfile& clone)
01773 : Profile<IMG_ELM_TYPE, MSK_ELM_TYPE>(clone.setting_),
01774 normFac(1.0/IMAGE_DATATYPE_MAXVAL),
01775 frameSeq_(NULL), maskSeq_(NULL),
01776 wSize_(0), width_(0), height_(0), nChl_(0), nElm_(0),
01777 pxls_(), maxElmIdx_(), neighbors_(),
01778 hue_(NULL), hueVariation_(NULL), hueVecX_(NULL), hueVecY_(NULL),
01779 lightness_(NULL), litVariation_(NULL), saturation_(NULL),
01780 satVariation_(NULL), intensity_(NULL), locationX_(NULL), locationY_(NULL),
01781 locationZ_(NULL), nPixels_(NULL), nElm_old_(0), nElm_front_(0)
01782
01783 {
01784 *this = clone;
01785 }
01786
01787 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
01788 void ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::cleanup()
01789 {
01790 if (hue_ != NULL) {
01791 delete [] hue_;
01792 hue_ = NULL;
01793 }
01794 if (hueVariation_ != NULL) {
01795 delete [] hueVariation_;
01796 hueVariation_ = NULL;
01797 }
01798 if (hueVecX_ != NULL) {
01799 delete [] hueVecX_;
01800 hueVecX_ = NULL;
01801 }
01802 if (hueVecY_ != NULL) {
01803 delete [] hueVecY_;
01804 hueVecY_ = NULL;
01805 }
01806 if (lightness_ != NULL) {
01807 delete [] lightness_;
01808 lightness_ = NULL;
01809 }
01810 if (litVariation_ != NULL) {
01811 delete [] litVariation_;
01812 litVariation_ = NULL;
01813 }
01814 if (saturation_ != NULL) {
01815 delete [] saturation_;
01816 saturation_ = NULL;
01817 }
01818 if (satVariation_ != NULL) {
01819 delete [] satVariation_;
01820 satVariation_ = NULL;
01821 }
01822 if (intensity_ != NULL) {
01823 for (unsigned int i=0; i < nElm_; i++) {
01824 if (intensity_[i] != NULL) delete [] intensity_[i];
01825 }
01826 delete [] intensity_;
01827 intensity_ = NULL;
01828 }
01829 if (locationX_ != NULL) {
01830 delete [] locationX_;
01831 locationX_ = NULL;
01832 }
01833 if (locationY_ != NULL) {
01834 delete [] locationY_;
01835 locationY_ = NULL;
01836 }
01837 if (locationZ_ != NULL) {
01838 delete [] locationZ_;
01839 locationZ_ = NULL;
01840 }
01841 if (nPixels_ != NULL) {
01842 delete [] nPixels_;
01843 nPixels_ = NULL;
01844 }
01845 }
01846
01847 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
01848 ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::~ElementProfile()
01849 {
01850 cleanup();
01851 }
01852
01853 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
01854 ElementProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>&
01855 ElementProfile<IMG_ELM_TYPE,MSK_ELM_TYPE>::
01856 operator=(const ElementProfile& clone)
01857 {
01858
01859 if (this == &clone) return *this;
01860
01861
01862 cleanup();
01863
01864
01865 this->similarity_ = clone.similarity_;
01866
01867
01868 this->indices_ = clone.indices_;
01869 this->identifiers_ = clone.identifiers_;
01870
01871
01872 this->frameSeq_ = clone.frameSeq_;
01873 this->maskSeq_ = clone.maskSeq_;
01874
01875 this->width_ = clone.width_;
01876 this->height_ = clone.height_;
01877 this->nChl_ = clone.nChl_;
01878 this->wSize_ = clone.wSize_;
01879 this->nElm_ = clone.nElm_;
01880
01881 this->pxls_ = clone.pxls_;
01882 this->maxElmIdx_ = clone.maxElmIdx_;
01883
01884 this->neighbors_ = clone.neighbors_;
01885
01886
01887 if (clone.hue_ != NULL) {
01888 this->hue_ = new double [clone.nElm_];
01889 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
01890 this->hue_[elm] = clone.hue_[elm];
01891 }
01892 }
01893 else {
01894 this->hue_ = NULL;
01895 }
01896
01897 if (clone.hueVariation_ != NULL) {
01898 this->hueVariation_ = new double [clone.nElm_];
01899 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
01900 this->hueVariation_[elm] = clone.hueVariation_[elm];
01901 }
01902 }
01903 else {
01904 this->hueVariation_ = NULL;
01905 }
01906
01907 if (clone.hueVecX_ != NULL) {
01908 this->hueVecX_ = new double [clone.nElm_];
01909 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
01910 this->hueVecX_[elm] = clone.hueVecX_[elm];
01911 }
01912 }
01913 else {
01914 this->hueVecX_ = NULL;
01915 }
01916
01917 if (clone.hueVecY_ != NULL) {
01918 this->hueVecY_ = new double [clone.nElm_];
01919 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
01920 this->hueVecY_[elm] = clone.hueVecY_[elm];
01921 }
01922 }
01923 else {
01924 this->hueVecY_ = NULL;
01925 }
01926
01927 if (clone.lightness_ != NULL) {
01928 this->lightness_ = new double [clone.nElm_];
01929 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
01930 this->lightness_[elm] = clone.lightness_[elm];
01931 }
01932 }
01933 else {
01934 this->lightness_ = NULL;
01935 }
01936
01937 if (clone.litVariation_ != NULL) {
01938 this->litVariation_ = new double [clone.nElm_];
01939 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
01940 this->litVariation_[elm] = clone.litVariation_[elm];
01941 }
01942 }
01943 else {
01944 this->litVariation_ = NULL;
01945 }
01946
01947 if (clone.saturation_ != NULL) {
01948 this->saturation_ = new double [clone.nElm_];
01949 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
01950 this->saturation_[elm] = clone.saturation_[elm];
01951 }
01952 }
01953 else {
01954 this->saturation_ = NULL;
01955 }
01956
01957 if (clone.satVariation_ != NULL) {
01958 this->satVariation_ = new double [clone.nElm_];
01959 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
01960 this->satVariation_[elm] = clone.satVariation_[elm];
01961 }
01962 }
01963 else {
01964 this->satVariation_ = NULL;
01965 }
01966
01967 if (clone.intensity_ != NULL) {
01968 this->intensity_ = new double* [clone.nElm_];
01969 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
01970 if (clone.intensity_[elm] == NULL) {
01971 this->intensity_[elm] = NULL;
01972 continue;
01973 }
01974 else {
01975 this->intensity_[elm] = new double [clone.frameSeq_.nChl()];
01976 for (unsigned int chl=0; chl<clone.frameSeq_.nChl(); chl++)
01977 this->intensity_[elm][chl] = clone.intensity_[elm][chl];
01978 }
01979 }
01980 }
01981 else {
01982 this->intensity_ = NULL;
01983 }
01984
01985 if (clone.locationX_ != NULL) {
01986 this->locationX_ = new double [clone.nElm_];
01987 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
01988 this->locationX_[elm] = clone.locationX_[elm];
01989 }
01990 }
01991 else {
01992 this->locationX_ = NULL;
01993 }
01994
01995 if (clone.locationY_ != NULL) {
01996 this->locationY_ = new double [clone.nElm_];
01997 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
01998 this->locationY_[elm] = clone.locationY_[elm];
01999 }
02000 }
02001 else {
02002 this->locationY_ = NULL;
02003 }
02004
02005 if (clone.locationZ_ != NULL) {
02006 this->locationZ_ = new double [clone.nElm_];
02007 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
02008 this->locationZ_[elm] = clone.locationZ_[elm];
02009 }
02010 }
02011 else {
02012 this->locationZ_ = NULL;
02013 }
02014
02015 if (clone.nPixels_ != NULL) {
02016 this->nPixels_ = new unsigned int [clone.nElm_];
02017 for (unsigned int elm=0; elm<clone.nElm_; elm++) {
02018 this->nPixels_[elm] = clone.nPixels_[elm];
02019 }
02020 }
02021 else {
02022 this->nPixels_ = NULL;
02023 }
02024
02025
02026 this->setting_ = clone.setting_;
02027
02028 return *this;
02029 }
02030
02031 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02032 const MSK_ELM_TYPE& SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
02033 mask(const MSK_ELM_TYPE& elmId) const
02034 {
02035 typename std::map<MSK_ELM_TYPE, MSK_ELM_TYPE>::const_iterator mask_p;
02036
02037 mask_p = mask_.find(elmId);
02038 if (mask_p == mask_.end())
02039 throw ncutException("element not found");
02040
02041 return mask_p->second;
02042 }
02043
02044 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02045 unsigned int SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
02046 age(const MSK_ELM_TYPE& segId) const
02047 {
02048 typename std::map<MSK_ELM_TYPE, unsigned int>::const_iterator age_p;
02049
02050 age_p = age_.find(segId);
02051 if (age_p == age_.end())
02052 throw ncutException("segment not found");
02053
02054 return age_p->second;
02055 }
02056
02057 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02058 bool SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
02059 inactive(const MSK_ELM_TYPE& segId) const
02060 {
02061 typename std::map<MSK_ELM_TYPE, bool>::const_iterator inactive_p;
02062
02063 inactive_p = inactive_.find(segId);
02064 if (inactive_p == inactive_.end())
02065 throw ncutException("segment not found");
02066
02067 return inactive_p->second;
02068 }
02069
02070 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02071 const MSK_ELM_TYPE& SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
02072 simSeg(const MSK_ELM_TYPE& segId) const
02073 {
02074 typename std::map<MSK_ELM_TYPE, MSK_ELM_TYPE>::const_iterator simSeg_p;
02075
02076 simSeg_p = simSeg_.find(segId);
02077 if (simSeg_p == simSeg_.end())
02078 throw ncutException("segment not found");
02079
02080 return simSeg_p->second;
02081 }
02082
02083 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02084 double SegmentProfile<IMG_ELM_TYPE,MSK_ELM_TYPE>::
02085 dif(unsigned int idx1, unsigned int idx2,
02086 const std::vector<std::vector<MSK_ELM_TYPE> >* elements) const
02087 {
02088 if (idx1 >= nSeg_ || idx2 >= nSeg_)
02089 throw ncutException("segment index out of range");
02090 if (idx1 == idx2)
02091 throw ncutException("segment indices must be different");
02092
02093 if (elements == NULL)
02094 elements = &this->elements_;
02095
02096
02097 bool stretchSat = (this->setting_->stretch_saturation_merge == 1);
02098
02099
02100 double point0_x = this->setting_->stretch_point0_x;
02101 double point0_y = this->setting_->stretch_point0_y;
02102 double point1_x = this->setting_->stretch_point1_x;
02103 double point1_y = this->setting_->stretch_point1_y;
02104 double point2_x = this->setting_->stretch_point2_x;
02105 double point2_y = this->setting_->stretch_point2_y;
02106 double point3_x = this->setting_->stretch_point3_x;
02107 double point3_y = this->setting_->stretch_point3_y;
02108
02109 if (point0_x < 0.0 || point0_x > IMAGE_DATATYPE_MAXVAL)
02110 throw ncutException("invalid stretch point 0 coordinate x");
02111 if (point0_y < 0.0 || point0_y > IMAGE_DATATYPE_MAXVAL)
02112 throw ncutException("invalid stretch point 0 coordinate y");
02113 if (point1_x <= point0_x || point1_x > IMAGE_DATATYPE_MAXVAL)
02114 throw ncutException("invalid stretch point 1 coordinate x");
02115 if (point1_y < 0.0 || point1_y > IMAGE_DATATYPE_MAXVAL)
02116 throw ncutException("invalid stretch point 1 coordinate y");
02117 if (point2_x <= point1_x || point2_x > IMAGE_DATATYPE_MAXVAL)
02118 throw ncutException("invalid stretch point 2 coordinate x");
02119 if (point2_y < 0.0 || point2_y > IMAGE_DATATYPE_MAXVAL)
02120 throw ncutException("invalid stretch point 2 coordinate y");
02121 if (point3_x <= point2_x || point3_x > IMAGE_DATATYPE_MAXVAL)
02122 throw ncutException("invalid stretch point 3 coordinate x");
02123 if (point3_y < 0.0 || point3_y > IMAGE_DATATYPE_MAXVAL)
02124 throw ncutException("invalid stretch point 3 coordinate y");
02125
02126
02127 const double perctFac = 100.0*elmProfile_->normFac;
02128
02129
02130 unsigned int neighborPxls = 0;
02131 typename std::vector<MSK_ELM_TYPE>::const_iterator elm1;
02132 typename std::vector<MSK_ELM_TYPE>::const_iterator elm2;
02133
02134 for (elm1 = (*elements)[idx1].begin();
02135 elm1 != (*elements)[idx1].end(); elm1++) {
02136 unsigned int elm1Idx = elmProfile_->idx(*elm1);
02137
02138 for (elm2 = (*elements)[idx2].begin();
02139 elm2 != (*elements)[idx2].end(); elm2++) {
02140 unsigned int elm2Idx = elmProfile_->idx(*elm2);
02141
02142
02143 if (elmProfile_->neighbors()->val(elm1Idx,elm2Idx) > 0.5)
02144 neighborPxls += static_cast<unsigned int>
02145 (elmProfile_->neighbors()->val(elm1Idx,elm2Idx));
02146 }
02147 }
02148
02149
02150 if (neighborPxls == 0)
02151 return 500.0;
02152
02153
02154
02155 double cutIntDif = 0.0;
02156 unsigned int pxlCount = 0;
02157
02158 for (elm1 = (*elements)[idx1].begin();
02159 elm1 != (*elements)[idx1].end(); elm1++) {
02160 unsigned int elm1Idx = elmProfile_->idx(*elm1);
02161
02162 for (elm2 = (*elements)[idx2].begin();
02163 elm2 != (*elements)[idx2].end(); elm2++) {
02164 unsigned int elm2Idx = elmProfile_->idx(*elm2);
02165
02166
02167 double hueVecX_i, hueVecY_i, hueVecX_j, hueVecY_j;
02168 if (stretchSat) {
02169 double sat_i = elmProfile_->saturation(elm1Idx);
02170 double sat_j = elmProfile_->saturation(elm2Idx);
02171
02172 double interPos;
02173 double stretch_i;
02174 if (sat_i <= point1_x) {
02175 interPos = (sat_i-point0_x) / (point1_x-point0_x);
02176 stretch_i = (interPos*point1_y) + (1-interPos)*point0_y;
02177 }
02178 else if (sat_i <= point2_x) {
02179 interPos = (sat_i-point1_x) / (point2_x-point1_x);
02180 stretch_i = (interPos*point2_y) + (1-interPos)*point1_y;
02181 }
02182 else {
02183 interPos = (sat_i-point2_x) / (point3_x-point2_x);
02184 stretch_i = (interPos*point3_y) + (1-interPos)*point2_y;
02185 }
02186
02187 double stretch_j;
02188 if (sat_j <= point1_x) {
02189 interPos = (sat_j-point0_x) / (point1_x-point0_x);
02190 stretch_j = (interPos*point1_y) + (1-interPos)*point0_y;
02191 }
02192 else if (sat_j <= point2_x) {
02193 interPos = (sat_j-point1_x) / (point2_x-point1_x);
02194 stretch_j = (interPos*point2_y) + (1-interPos)*point1_y;
02195 }
02196 else {
02197 interPos = (sat_j-point2_x) / (point3_x-point2_x);
02198 stretch_j = (interPos*point3_y) + (1-interPos)*point2_y;
02199 }
02200
02201 if (sat_i > 0.0) {
02202 hueVecX_i = elmProfile_->hueVecX(elm1Idx) * (stretch_i/sat_i);
02203 hueVecY_i = elmProfile_->hueVecY(elm1Idx) * (stretch_i/sat_i);
02204 }
02205 else {
02206 hueVecX_i = 0.0;
02207 hueVecY_i = 0.0;
02208 }
02209 if (sat_j > 0.0) {
02210 hueVecX_j = elmProfile_->hueVecX(elm2Idx) * (stretch_j/sat_j);
02211 hueVecY_j = elmProfile_->hueVecY(elm2Idx) * (stretch_j/sat_j);
02212 }
02213 else {
02214 hueVecX_j = 0.0;
02215 hueVecY_j = 0.0;
02216 }
02217 }
02218 else {
02219 hueVecX_i = elmProfile_->hueVecX(elm1Idx);
02220 hueVecY_i = elmProfile_->hueVecY(elm1Idx);
02221 hueVecX_j = elmProfile_->hueVecX(elm2Idx);
02222 hueVecY_j = elmProfile_->hueVecY(elm2Idx);
02223 }
02224
02225
02226 double hueSatDifX = fabs(hueVecX_i - hueVecX_j);
02227 hueSatDifX *= perctFac;
02228
02229 double hueSatDifY = fabs(hueVecY_i - hueVecY_j);
02230 hueSatDifY *= perctFac;
02231
02232 double litDif = fabs(elmProfile_->lightness(elm1Idx) -
02233 elmProfile_->lightness(elm2Idx));
02234 litDif *= perctFac;
02235
02236
02237 double intDif = sqrt(hueSatDifX*hueSatDifX +
02238 hueSatDifY*hueSatDifY +
02239 litDif*litDif);
02240
02241
02242 unsigned int pxls = elmProfile_->nPixels(elm1Idx) +
02243 elmProfile_->nPixels(elm2Idx);
02244 cutIntDif += intDif * pxls;
02245 pxlCount += pxls;
02246 }
02247 }
02248
02249 return cutIntDif / pxlCount;
02250 }
02251
02252 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02253 int SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::calcSim()
02254 {
02255 this->similarity_.changeDimension(nSeg_);
02256
02257 try {
02258 for (unsigned int seg1=0; seg1<nSeg_; seg1++) {
02259 for (unsigned int seg2=seg1+1; seg2<nSeg_; seg2++) {
02260 this->similarity_.val(seg2,seg1) = this->similarity_.val(seg1,seg2) =
02261 500.0 - dif(seg1, seg2);
02262 }
02263 }
02264 }
02265 catch (ncutException& e) {
02266 this->similarity_.changeDimension(0);
02267 throw;
02268 }
02269
02270 return 0;
02271 }
02272
02273 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02274 void SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
02275 attachIds(const SegmentProfile<IMG_ELM_TYPE,MSK_ELM_TYPE>& last)
02276 {
02277 unsigned int adultAge = this->setting_->segment_adult_age;
02278 double ageWeight = this->setting_->segment_age_weight;
02279 double attachSim = this->setting_->min_segment_attach_similarity;
02280
02281 if (adultAge == 0)
02282 throw ncutException("invalid segment adult age, must be >= 1.");
02283 if (ageWeight < 0.0 || ageWeight > 1.0)
02284 throw ncutException("invalid segment age weight, must be in [0,1].");
02285 if (attachSim < 0.0 || attachSim > 1.0)
02286 throw ncutException(
02287 "invalid minimum segment attach similarity, must be in [0,1].");
02288
02289
02290
02291 double similarity[last.elements_.size()][elements_.size()];
02292 unsigned int totalPxlArea_last[last.elements_.size()];
02293 unsigned int totalPxlArea[elements_.size()];
02294
02295 for (unsigned int seg_last=0; seg_last<last.elements_.size(); seg_last++)
02296 for (unsigned int seg=0; seg<elements_.size(); seg++)
02297 similarity[seg_last][seg] = 0.0;
02298
02299 for (unsigned int seg_last=0; seg_last<last.elements_.size(); seg_last++)
02300 totalPxlArea_last[seg_last] = 0;
02301
02302 for (unsigned int seg=0; seg<elements_.size(); seg++)
02303 totalPxlArea[seg] = 0;
02304
02305
02306
02307 typename std::vector<std::vector<MSK_ELM_TYPE> >::iterator seg;
02308 typename std::vector<MSK_ELM_TYPE>::iterator elmId;
02309 unsigned int segIdx;
02310 unsigned long pxlArea;
02311 for (seg = elements_.begin(), segIdx=0;
02312 seg != elements_.end(); seg++, segIdx++) {
02313 for(elmId = seg->begin(); elmId != seg->end(); elmId++) {
02314 unsigned int elmIdx = elmProfile_->idx(*elmId);
02315
02316 if (static_cast<unsigned int>(elmProfile_->locationZ(elmIdx)) <
02317 elmProfile_->wSize()-1) {
02318 unsigned int segIdx_last = last.idx(last.mask(*elmId));
02319
02320 pxlArea = elmProfile_->nPixels(elmIdx);
02321 similarity[segIdx_last][segIdx] += pxlArea;
02322 totalPxlArea_last[segIdx_last] += pxlArea;
02323 totalPxlArea[segIdx] += pxlArea;
02324 }
02325 }
02326 }
02327
02328
02329
02330 double age;
02331 for (unsigned int seg_last=0; seg_last<last.elements_.size(); seg_last++) {
02332 for (unsigned int seg=0; seg<elements_.size(); seg++) {
02333
02334
02335
02336
02337
02338
02339 similarity[seg_last][seg] *= 2.0;
02340 similarity[seg_last][seg] /= (totalPxlArea_last[seg_last] +
02341 totalPxlArea[seg]);
02342
02343
02344
02345
02346
02347
02348 age = last.age(last.id(seg_last));
02349 age /= adultAge;
02350 if (age > 1.0) age = 1.0;
02351 age = (1.0-ageWeight)+(age*ageWeight);
02352 similarity[seg_last][seg] *= age;
02353 }
02354 }
02355
02356
02357
02358
02359
02360
02361
02362
02363
02364
02365
02366
02367
02368
02369
02370
02371
02372
02373 this->identifiers_.resize(nSeg_);
02374 unsigned int similarSegment;
02375 std::vector<double> maxSim(nSeg_,0.0);
02376 for (unsigned int seg=0; seg<nSeg_; seg++) {
02377
02378 double currentMaxSim = -1.0;
02379 for (unsigned int seg_last=0; seg_last<last.nSeg_; seg_last++) {
02380 if (similarity[seg_last][seg] > currentMaxSim) {
02381 currentMaxSim = similarity[seg_last][seg];
02382 similarSegment = seg_last;
02383 }
02384 }
02385
02386
02387
02388
02389
02390
02391
02392 if (currentMaxSim > attachSim) {
02393
02394 MSK_ELM_TYPE id = last.id(similarSegment);
02395 typename std::map<MSK_ELM_TYPE, unsigned int>::iterator otherSeg_p =
02396 this->indices_.find(id);
02397
02398 if (otherSeg_p == this->indices_.end()) {
02399
02400
02401
02402 this->identifiers_[seg] = id;
02403 this->indices_[id] = seg;
02404 maxSim[seg] = currentMaxSim;
02405 continue;
02406 }
02407 else {
02408
02409 int otherSeg = otherSeg_p->second;
02410
02411
02412 if (currentMaxSim > maxSim[otherSeg]) {
02413
02414
02415
02416 this->identifiers_[seg] = id;
02417 this->indices_[id] = seg;
02418 maxSim[seg] = currentMaxSim;
02419
02420
02421 if (nextSegId_ > MAX_SEG_ID)
02422 throw ncutCritical("not enough segment Id numbers");
02423 this->identifiers_[otherSeg] = nextSegId_;
02424 this->indices_[nextSegId_] = otherSeg;
02425 nextSegId_++;
02426 continue;
02427 }
02428
02429
02430 }
02431 }
02432
02433
02434 if (nextSegId_ > MAX_SEG_ID)
02435 throw ncutCritical("not enough segment Id numbers");
02436
02437 this->identifiers_[seg] = nextSegId_;
02438 this->indices_[nextSegId_] = seg;
02439 nextSegId_++;
02440
02441
02442
02443
02444
02445 }
02446 }
02447
02448 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02449 void SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
02450 attach(const SegmentProfile<IMG_ELM_TYPE,MSK_ELM_TYPE>& last)
02451 {
02452 unsigned int segDelay = this->setting_->segment_delay;
02453 unsigned int adultAge = this->setting_->segment_adult_age;
02454
02455 if (segDelay < 0)
02456 throw ncutException("invalid segment delay, must be >= 0.");
02457 if (adultAge == 0)
02458 throw ncutException("invalid segment adult age, must be >= 1.");
02459
02460 if (last.nSeg_ == 0) {
02461
02462
02463
02464 typename std::vector<MSK_ELM_TYPE>::iterator elm;
02465 for (; nextSegId_<elements_.size(); nextSegId_++) {
02466 this->identifiers_.push_back(static_cast<MSK_ELM_TYPE>(nextSegId_));
02467 this->indices_[this->id(nextSegId_)] = nextSegId_;
02468 for (elm = elements_[nextSegId_].begin();
02469 elm != elements_[nextSegId_].end(); elm++) {
02470 mask_[*elm] = this->id(nextSegId_);
02471 }
02472 age_[this->id(nextSegId_)] = (adultAge>segDelay) ? adultAge+1 : segDelay+1;
02473 inactive_[this->id(nextSegId_)] = false;
02474 simSeg_[this->id(nextSegId_)] = this->id(nextSegId_);
02475 }
02476 }
02477 else {
02478
02479
02480
02481 attachIds(last);
02482
02483
02484
02485
02486
02487
02488
02489 typename std::vector<MSK_ELM_TYPE>::iterator elm;
02490 for (unsigned int seg=0; seg<nSeg_; seg++) {
02491 for (elm = elements_[seg].begin(); elm != elements_[seg].end(); elm++)
02492 mask_[*elm] = this->id(seg);
02493 }
02494
02495
02496 typename std::map<MSK_ELM_TYPE, bool>::iterator inactive;
02497 for (inactive = inactive_.begin(); inactive!=inactive_.end(); inactive++)
02498 inactive->second = true;
02499 typename std::map<MSK_ELM_TYPE, unsigned int>::const_iterator last_active;
02500 for (unsigned int seg=0; seg < nSeg_; seg++) {
02501 MSK_ELM_TYPE segId = this->id(seg);
02502 last_active = last.age_.find(segId);
02503 if (last_active == last.age_.end()) {
02504 age_[segId] = 1;
02505 inactive_[segId] = false;
02506 }
02507 else {
02508 age_[segId] += last_active->second + 1;
02509 inactive_[segId] = false;
02510 }
02511
02512 }
02513
02514
02515
02516
02517
02518
02519
02520
02521
02522
02523
02524
02525
02526
02527
02528
02529 for (unsigned int seg=0; seg<nSeg_; seg++) {
02530 if (age_[this->id(seg)] < segDelay) {
02531
02532
02533
02534 double minDif = 600.0;
02535 int minDifSeg = -1;
02536 for (unsigned int seg2=0; seg2<nSeg_; seg2++) {
02537 if (seg == seg2) continue;
02538 if (age_[this->id(seg2)] < segDelay) continue;
02539 double difference = dif(seg, seg2);
02540 if (difference < minDif) {
02541 minDifSeg = seg2;
02542 minDif = difference;
02543 }
02544 }
02545 if (minDifSeg == -1)
02546 simSeg_[this->id(seg)] = this->id(seg);
02547 else
02548 simSeg_[this->id(seg)] = this->id(minDifSeg);
02549 }
02550 else {
02551 simSeg_[this->id(seg)] = this->id(seg);
02552 }
02553 }
02554
02555 }
02556 }
02557
02558 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02559 int SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
02560 push(const ElementProfile<IMG_ELM_TYPE,MSK_ELM_TYPE>* elmProfilee,
02561 const std::vector<std::vector<MSK_ELM_TYPE> >& segmentss)
02562 {
02563 if (elmProfilee == NULL)
02564 throw ncutException("sp invalid elment profile (NULL)");
02565 if (elmProfilee->nElm() == 0)
02566 throw ncutException("invalid elment profile (does not contain elements)");
02567 if (segmentss.empty())
02568 throw ncutException("invalid segments list (empty)");
02569
02570
02571 SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE> last(*this);
02572
02573
02574 cleanup();
02575 this->indices_.clear();
02576 this->identifiers_.clear();
02577 mask_.clear();
02578 idxMask_.clear();
02579 age_.clear();
02580 simSeg_.clear();
02581
02582
02583 elmProfile_ = elmProfilee;
02584 nChl_ = elmProfilee->nChl();
02585 nSeg_ = segmentss.size();
02586 elements_ = segmentss;
02587 idxMask_.resize(elmProfilee->nElm());
02588 neighbors_.changeDimension(nSeg_);
02589 intensity_ = new double* [nSeg_];
02590 for (unsigned int seg=0; seg<nSeg_; seg++)
02591 intensity_[seg] = new double[nChl_];
02592 locationX_ = new double [nSeg_];
02593 locationY_ = new double [nSeg_];
02594 nElements_ = new unsigned int [nSeg_];
02595 nPixels_ = new unsigned int [nSeg_];
02596
02597
02598 for (unsigned int seg=0; seg<nSeg_; seg++) {
02599 for (unsigned int chl = 0; chl<nChl_; chl++) {
02600 intensity_[seg][chl] = 0.0;
02601 }
02602 locationX_[seg] = 0.0;
02603 locationY_[seg] = 0.0;
02604 nElements_[seg] = 0;
02605 nPixels_[seg] = 0;
02606 }
02607
02608
02609 try {
02610 unsigned int elmIdx;
02611 typename std::vector<MSK_ELM_TYPE>::iterator elmId_p;
02612 for (unsigned int segIdx=0; segIdx<elements_.size(); segIdx++) {
02613 for (elmId_p = elements_[segIdx].begin();
02614 elmId_p != elements_[segIdx].end(); elmId_p++) {
02615
02616
02617 elmIdx = elmProfile_->idx(*elmId_p);
02618
02619
02620 idxMask_[elmIdx] = segIdx;
02621
02622
02623 for (unsigned int chl = 0; chl<nChl_; chl++)
02624 intensity_[segIdx][chl] +=
02625 elmProfile_->intensity(elmIdx,chl) *
02626 elmProfile_->nPixels(elmIdx);
02627 nElements_[segIdx]++;
02628 nPixels_[segIdx] += elmProfile_->nPixels(elmIdx);
02629
02630
02631 locationX_[segIdx] += elmProfile_->locationX(elmIdx) *
02632 elmProfile_->nPixels(elmIdx);
02633 locationY_[segIdx] += elmProfile_->locationY(elmIdx) *
02634 elmProfile_->nPixels(elmIdx);
02635 }
02636 }
02637 }
02638 catch (ncutException& e) {
02639 *this = last;
02640 throw;
02641 }
02642
02643
02644
02645
02646 for (unsigned int seg=0; seg<nSeg_; seg++) {
02647 if (nPixels_[seg] == 0) {
02648 for (unsigned int chl = 0; chl<nChl_; chl++)
02649 intensity_[seg][chl] = 0.0;
02650 locationX_[seg] = -1;
02651 locationY_[seg] = -1;
02652
02653 continue;
02654 }
02655
02656 for (unsigned int chl=0; chl<nChl_; chl++) {
02657 intensity_[seg][chl] /= nPixels_[seg];
02658 }
02659
02660 locationX_[seg] /= nPixels_[seg];
02661 locationY_[seg] /= nPixels_[seg];
02662 }
02663
02664
02665 for (unsigned int i=0; i<elmProfile_->nElm(); i++) {
02666 for (unsigned int j=i+1; j<elmProfile_->nElm(); j++) {
02667 if (idxMask_[i] != idxMask_[j] &&
02668 elmProfile_->neighbors()->val(i,j) > 0.5) {
02669 neighbors_.val(idxMask_[i],idxMask_[j]) +=
02670 elmProfile_->neighbors()->val(i,j);
02671 neighbors_.val(idxMask_[j],idxMask_[i]) =
02672 neighbors_.val(idxMask_[i],idxMask_[j]);
02673 }
02674 }
02675 }
02676
02677 try {
02678 attach(last);
02679 }
02680 catch (ncutException& e) {
02681 *this = last;
02682 throw;
02683 }
02684
02685 return 0;
02686 }
02687
02688 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02689 SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
02690 SegmentProfile(const SegmentProfile& clone)
02691 : Profile<IMG_ELM_TYPE, MSK_ELM_TYPE>(clone.setting_),
02692 elements_(),
02693 elmProfile_(NULL),
02694 mask_(),
02695 idxMask_(),
02696 age_(),
02697 inactive_(),
02698 simSeg_(),
02699 nSeg_(0),
02700 nChl_(0),
02701 nextSegId_(0),
02702 neighbors_(),
02703 intensity_(NULL),
02704 locationX_(NULL),
02705 locationY_(NULL),
02706 nElements_(NULL),
02707 nPixels_(NULL)
02708 {
02709 *this = clone;
02710 }
02711
02712 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02713 SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
02714 SegmentProfile(const Setting* setting)
02715 : Profile<IMG_ELM_TYPE, MSK_ELM_TYPE>(setting),
02716 elements_(),
02717 elmProfile_(NULL),
02718 mask_(),
02719 idxMask_(),
02720 age_(),
02721 inactive_(),
02722 simSeg_(),
02723 nSeg_(0),
02724 nChl_(0),
02725 nextSegId_(0),
02726 neighbors_(),
02727 intensity_(NULL),
02728 locationX_(NULL),
02729 locationY_(NULL),
02730 nElements_(NULL),
02731 nPixels_(NULL)
02732 {
02733 return;
02734 }
02735
02736 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02737 void SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::cleanup()
02738 {
02739 if (nElements_ != NULL) {
02740 delete [] nElements_;
02741 nElements_ = NULL;
02742 }
02743 if (nPixels_ != NULL) {
02744 delete [] nPixels_;
02745 nPixels_ = NULL;
02746 }
02747 if (intensity_ != NULL) {
02748 for (unsigned int seg=0; seg < nSeg_; seg++) {
02749 if (intensity_[seg] != NULL) delete [] intensity_[seg];
02750 }
02751 delete [] intensity_;
02752 intensity_ = NULL;
02753 }
02754 if (locationX_ != NULL) {
02755 delete [] locationX_;
02756 locationX_ = NULL;
02757 }
02758 if (locationY_ != NULL) {
02759 delete [] locationY_;
02760 locationY_ = NULL;
02761 }
02762 }
02763
02764 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02765 SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::~SegmentProfile()
02766 {
02767 cleanup();
02768 }
02769
02770 template <class IMG_ELM_TYPE, class MSK_ELM_TYPE>
02771 SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>&
02772 SegmentProfile<IMG_ELM_TYPE, MSK_ELM_TYPE>::
02773 operator=(const SegmentProfile& clone)
02774 {
02775
02776 if (this == &clone) return *this;
02777
02778
02779 cleanup();
02780
02781
02782 this->indices_ = clone.indices_;
02783 this->identifiers_ = clone.identifiers_;
02784 this->similarity_ = clone.similarity_;
02785 this->elements_ = clone.elements_;
02786 this->elmProfile_ = clone.elmProfile_;
02787 this->mask_ = clone.mask_;
02788 this->idxMask_ = clone.idxMask_;
02789 this->age_ = clone.age_;
02790 this->inactive_ = clone.inactive_;
02791 this->simSeg_ = clone.simSeg_;
02792 this->nChl_ = clone.nChl_;
02793 this->nSeg_ = clone.nSeg_;
02794 this->nextSegId_ = clone.nextSegId_;
02795 this->neighbors_ = clone.neighbors_;
02796
02797
02798 if (clone.intensity_ != NULL) {
02799 this->intensity_ = new double* [clone.nSeg_];
02800 for (unsigned int seg=0; seg<clone.nSeg_; seg++) {
02801 if (clone.intensity_[seg] == NULL) {
02802 this->intensity_[seg] = NULL;
02803 continue;
02804 }
02805 else {
02806 this->intensity_[seg] = new double [clone.nChl_];
02807 for (unsigned int chl=0; chl<clone.nChl_; chl++) {
02808 this->intensity_[seg][chl] = clone.intensity_[seg][chl];
02809 }
02810 }
02811 }
02812 }
02813 else {
02814 this->intensity_ = NULL;
02815 }
02816
02817 if (clone.locationX_ != NULL) {
02818 this->locationX_ = new double [clone.nSeg_];
02819 for (unsigned int seg=0; seg<clone.nSeg_; seg++) {
02820 this->locationX_[seg] = clone.locationX_[seg];
02821 }
02822 }
02823 else {
02824 this->locationX_ = NULL;
02825 }
02826
02827 if (clone.locationY_ != NULL) {
02828 this->locationY_ = new double [clone.nSeg_];
02829 for (unsigned int seg=0; seg<clone.nSeg_; seg++) {
02830 this->locationY_[seg] = clone.locationY_[seg];
02831 }
02832 }
02833 else {
02834 this->locationY_ = NULL;
02835 }
02836
02837 if (clone.nElements_ != NULL) {
02838 this->nElements_ = new unsigned int [clone.nSeg_];
02839 for (unsigned int seg=0; seg<clone.nSeg_; seg++) {
02840 this->nElements_[seg] = clone.nElements_[seg];
02841 }
02842 }
02843 else {
02844 this->nElements_ = NULL;
02845 }
02846
02847 if (clone.nPixels_ != NULL) {
02848 this->nPixels_ = new unsigned int [clone.nSeg_];
02849 for (unsigned int seg=0; seg<clone.nSeg_; seg++) {
02850 this->nPixels_[seg] = clone.nPixels_[seg];
02851 }
02852 }
02853 else {
02854 this->nPixels_ = NULL;
02855 }
02856
02857
02858 this->setting_ = clone.setting_;
02859
02860 return *this;
02861 }
02862
02863 }
02864
02865 #endif // MOTIONPROFILE_H