Changeset 90


Ignore:
Timestamp:
07/22/2009 10:59:15 AM (10 years ago)
Author:
dave
Message:

PCA updated - projection and resynthesis seem to work

Location:
foam/trunk/vision
Files:
11 edited

Legend:

Unmodified
Added
Removed
  • foam/trunk/vision/Makefile

    r89 r90  
    88        src/LDAClassifier.cpp\ 
    99        src/PCA.cpp\ 
    10         src/SVD.cpp 
     10        src/SVD.cpp\ 
     11        src/tinyxml.cpp\ 
     12        src/tinyxmlerror.cpp\ 
     13        src/tinyxmlparser.cpp 
    1114 
    1215CCFLAGS = `pkg-config --cflags opencv` -ggdb -Wall -O3 -ffast-math -Wno-unused -DTIXML_USE_STL  
  • foam/trunk/vision/src/Classifier.h

    r86 r90  
    4040        typedef std::vector<Vector<float> > FeatureVec; 
    4141        typedef std::map<int,FeatureVec> FeatureMap; 
     42        typedef std::map<int,FeatureVec> GroupMap; 
    4243         
    4344        FeatureMap m_Features; 
     45        GroupMap m_Groups; 
    4446        unsigned int m_FeatureSize; 
    4547 
  • foam/trunk/vision/src/Image.cpp

    r89 r90  
    132132         
    133133        return ret; 
     134} 
     135 
     136Image &Image::operator=(const Image &other) 
     137{ 
     138        m_Image=cvCloneImage(other.m_Image); 
     139        return *this; 
    134140} 
    135141 
  • foam/trunk/vision/src/Image.h

    r89 r90  
    3636        Image operator-(const Image &other); 
    3737        Image operator+(const Image &other); 
     38        Image &operator=(const Image &other); 
    3839 
    3940        void PrintInfo(); 
  • foam/trunk/vision/src/LDAClassifier.cpp

    r86 r90  
    3333void LDAClassifier::CalcGroupMeans() 
    3434{        
    35         for (FeatureMap::iterator i=m_Features.begin(); 
    36                 i!=m_Features.end(); ++i) 
     35        for (GroupMap::iterator i=m_Groups.begin(); 
     36                i!=m_Groups.end(); ++i) 
    3737        { 
    3838                m_GroupMean[i->first]=Vector<float>(m_FeatureSize); 
     
    4747 
    4848void LDAClassifier::CalcMeanCorrected() 
    49 { 
     49{        
     50/*      CalcMean(); 
     51         
     52        // copy the training data :/ 
     53        m_MeanCorrectedGroups = m_Groups; 
     54         
     55        for (GroupMap::iterator i=m_MeanCorrectedGroups.begin(); 
     56                i!=m_MeanCorrectedGroups.end(); ++i) 
     57        { 
     58                Matrix<float> Group(i->second.size(), m_FeatureSize); 
     59                unsigned int count=0; 
     60                for (FeatureVec::iterator vi = i->second.begin(); vi!=i->second.end(); ++vi) 
     61                { 
     62                        Group.SetRowVector(count++,*vi); 
     63                } 
     64                m_MeanCorrectedGroups[i->first]=Group; 
     65        }*/ 
    5066} 
    5167 
    52 void LDAClassifier::CalcCovariance() 
     68void LDAClassifier::CalcGroupCovariance() 
    5369{ 
     70         
    5471} 
    5572 
  • foam/trunk/vision/src/LDAClassifier.h

    r86 r90  
    3939        void CalcGroupMeans(); 
    4040        void CalcMeanCorrected(); 
    41         void CalcCovariance(); 
     41        void CalcGroupCovariance(); 
    4242        void CalcPooledCovariance(); 
    4343        void CalcPriorProbablity(); 
    4444         
    4545        std::map<int,Vector<float> > m_GroupMean; 
    46         std::map<int,Matrix<float> > m_MeanCorrected; 
    47         std::map<int,Matrix<float> > m_Covariance; 
     46        std::map<int,Matrix<float> > m_MeanCorrectedGroups; 
     47        std::map<int,Matrix<float> > m_GroupCovariance; 
     48 
    4849        //Matrix<T> m_PooledCovariance; 
    4950        //Vector<T> m_PriorProbability; 
  • foam/trunk/vision/src/Matrix.h

    r89 r90  
    9999        void Zero() { SetAll(0); } 
    100100        bool IsInf(); 
    101         Matrix Transposed(); 
     101        Matrix Transposed() const; 
    102102 
    103103        Matrix &operator=(const Matrix &other); 
     
    105105        Matrix operator-(const Matrix &other) const; 
    106106        Matrix operator*(const Matrix &other) const; 
     107        Vector<T> operator*(const Vector<T> &other) const; 
     108        Vector<T> VecMulTransposed(const Vector<T> &other) const; 
    107109        Matrix &operator+=(const Matrix &other); 
    108110        Matrix &operator-=(const Matrix &other); 
    109111        Matrix &operator*=(const Matrix &other); 
    110          
     112 
    111113        void SortRows(Vector<T> &v); 
    112114        void SortCols(Vector<T> &v); 
     115         
     116        Matrix CropRows(unsigned int s, unsigned int e); 
     117        Matrix CropCols(unsigned int s, unsigned int e); 
     118 
     119        void Save(FILE *f); 
     120        void Load(FILE *f); 
    113121 
    114122        static void RunTests(); 
     
    218226 
    219227template<class T> 
    220 Matrix<T> Matrix<T>::Transposed() 
     228Matrix<T> Matrix<T>::Transposed() const 
    221229{ 
    222230        Matrix<T> copy(*this); 
     
    281289                        { 
    282290                                ret[i][j]+=(*this)[i][k]*other[k][j]; 
     291                        } 
     292                } 
     293        } 
     294        return ret; 
     295} 
     296 
     297template<class T> 
     298Vector<T> Matrix<T>::operator*(const Vector<T> &other) const 
     299{ 
     300        assert(m_Cols==other.Size()); 
     301         
     302        Vector<T> ret(m_Rows); 
     303         
     304        for (unsigned int i=0; i<m_Rows; i++) 
     305        { 
     306                for (unsigned int j=0; j<other.Size(); j++) 
     307                { 
     308                        ret[i]=0; 
     309                        for (unsigned int k=0; k<m_Cols; k++) 
     310                        { 
     311                                ret[i]+=(*this)[i][k]*other[k]; 
     312                        } 
     313                } 
     314        } 
     315        return ret; 
     316} 
     317 
     318template<class T> 
     319Vector<T> Matrix<T>::VecMulTransposed(const Vector<T> &other) const 
     320{ 
     321        assert(m_Rows==other.Size()); 
     322         
     323        Vector<T> ret(m_Cols); 
     324         
     325        for (unsigned int i=0; i<m_Cols; i++) 
     326        { 
     327                for (unsigned int j=0; j<other.Size(); j++) 
     328                { 
     329                        ret[i]=0; 
     330                        for (unsigned int k=0; k<m_Rows; k++) 
     331                        { 
     332                                ret[i]+=(*this)[k][i]*other[k]; 
    283333                        } 
    284334                } 
     
    411461        } 
    412462} 
    413  
     463         
     464template<class T> 
     465Matrix<T> Matrix<T>::CropRows(unsigned int s, unsigned int e) 
     466{ 
     467        assert(s<e); 
     468        assert(s<m_Rows); 
     469        assert(e<=m_Rows); 
     470         
     471        Matrix r(e-s,m_Cols); 
     472        unsigned int c=0; 
     473        for(unsigned int i=s; i<e; i++) 
     474        { 
     475                r.SetRowVector(c,GetRowVector(i)); 
     476                c++; 
     477        } 
     478         
     479        return r; 
     480} 
     481 
     482template<class T> 
     483Matrix<T> Matrix<T>::CropCols(unsigned int s, unsigned int e) 
     484{ 
     485        assert(s<e); 
     486        assert(s<m_Cols); 
     487        assert(e<=m_Cols); 
     488         
     489        Matrix r(m_Rows,e-s); 
     490        unsigned int c=0; 
     491        for(unsigned int i=s; i<e; i++) 
     492        { 
     493                r.SetColVector(c,GetColVector(i)); 
     494                c++; 
     495        } 
     496         
     497        return r; 
     498} 
     499 
     500template<class T> 
     501void Matrix<T>::Save(FILE* f) 
     502{ 
     503        int version = 1;         
     504        fwrite(&version,1,sizeof(version),f); 
     505        fwrite(&m_Rows,1,sizeof(m_Rows),f); 
     506        fwrite(&m_Cols,1,sizeof(m_Cols),f); 
     507        fwrite(m_Data,1,sizeof(T)*m_Rows*m_Cols,f); 
     508} 
     509 
     510template<class T> 
     511void Matrix<T>::Load(FILE* f) 
     512{ 
     513        int version;     
     514        fread(&version,sizeof(version),1,f); 
     515        fread(&m_Rows,sizeof(m_Rows),1,f); 
     516        fread(&m_Cols,sizeof(m_Cols),1,f); 
     517        m_Data=new T[m_Rows*m_Cols]; 
     518        fread(m_Data,sizeof(T)*m_Rows*m_Cols,1,f); 
     519} 
    414520 
    415521template<class T> 
    416522void Matrix<T>::RunTests() 
    417523{ 
     524        Vector<T>::RunTests(); 
     525 
    418526        Matrix<T> m(10,10); 
    419527        m.SetAll(0); 
     
    432540        Matrix<T> c=a*b; 
    433541        assert(c[0][0]==11 && c[1][0]==29); 
     542         
     543        // test matrix * vector 
     544        Vector<T> d(3); 
     545        d[0]=3; 
     546        d[1]=1; 
     547        d[2]=2; 
     548        Vector<T> e=a*d; 
     549        assert(e[0]==11 && e[1]==29); 
     550         
     551        Matrix<T> f=a.CropCols(1,3); 
     552        assert(f.GetRows()==2 && f.GetCols()==2 && f[0][0]==2); 
     553        Matrix<T> g=a.CropRows(0,1); 
     554        assert(g.GetRows()==1 && g.GetCols()==3 && g[0][0]==1); 
     555 
    434556} 
    435557 
  • foam/trunk/vision/src/PCA.cpp

    r89 r90  
    3535{ 
    3636        // calculate the mean 
    37         Vector<float> Mean(m_FeatureSize); 
    38         Mean.Zero(); 
     37        m_Mean.Zero(); 
    3938        for (FeatureVec::iterator vi = m_Features.begin(); vi!=m_Features.end(); ++vi) 
    4039        { 
    41                 Mean+=*vi; 
     40                m_Mean+=*vi; 
    4241        } 
    4342         
    44         Mean/=m_Features.size(); 
     43        m_Mean/=m_Features.size(); 
    4544                 
    4645        // subtract the mean 
     
    4847        for (FeatureVec::iterator vi = m_Features.begin(); vi!=m_Features.end(); ++vi) 
    4948        { 
    50                 SubMean.push_back(*vi-Mean); 
     49                SubMean.push_back(*vi-m_Mean); 
    5150        } 
    5251         
     
    7069        } 
    7170        m_EigenValues = SVD(m_EigenTransform); 
     71        m_EigenTransform=m_EigenTransform.Transposed(); 
     72} 
     73 
     74void PCA::Compress(unsigned int s, unsigned int e) 
     75{ 
     76        m_EigenTransform=m_EigenTransform.CropRows(s,e); 
     77} 
     78 
     79Vector<float> PCA::Project(Vector<float> v) 
     80{ 
     81        return m_EigenTransform*v; 
     82} 
     83         
     84Vector<float> PCA::Synth(Vector<float> v) 
     85{ 
     86        return m_Mean+m_EigenTransform.VecMulTransposed(v); 
     87} 
     88 
     89void PCA::Save(FILE *f) 
     90{ 
     91        int version = 2; 
     92        fwrite(&version,sizeof(version),1,f); 
     93        m_EigenTransform.Save(f); 
     94        m_EigenValues.Save(f); 
     95        m_Mean.Save(f); 
     96        fwrite(&m_FeatureSize,sizeof(m_FeatureSize),1,f); 
     97        for (unsigned int i=0; i<m_Features.size(); i++) 
     98        { 
     99                m_Features[i].Save(f); 
     100        } 
     101         
     102} 
     103 
     104void PCA::Load(FILE *f) 
     105{ 
     106        int version;     
     107        fread(&version,sizeof(version),1,f); 
     108        m_EigenTransform.Load(f); 
     109        m_EigenValues.Load(f); 
     110         
     111        if (version == 1) 
     112        { 
     113                m_EigenTransform=m_EigenTransform.Transposed(); 
     114        } 
     115         
     116        if (version>1)  
     117        { 
     118                m_Mean.Load(f);  
     119                fread(&m_FeatureSize,sizeof(m_FeatureSize),1,f); 
     120                for (unsigned int i=0; i<m_Features.size(); i++) 
     121                { 
     122                        m_Features[i].Load(f); 
     123                } 
     124        } 
    72125} 
    73126 
    74127void PCA::RunTests() 
    75128{ 
     129        Matrix<float>::RunTests(); 
    76130        PCA pca(2); 
    77131         
     
    99153         
    100154        pca.Calculate(); 
     155                 
     156        in[0]=.69; in[1]=.49; 
     157        Vector<float> out = pca.Project(in); 
     158        assert(feq(out[0],-0.82797f) && feq(out[1],-0.175115f)); 
     159         
    101160} 
    102161 
  • foam/trunk/vision/src/PCA.h

    r89 r90  
    3737        void AddFeature(Vector<float> v) { m_Features.push_back(v); } 
    3838        void Calculate(); 
     39         
     40        // remove eigenvectors from the transform 
     41        void Compress(unsigned int s, unsigned int e); 
     42         
     43        Vector<float> Project(Vector<float> v); 
     44        Vector<float> Synth(Vector<float> v); 
    3945 
    4046        static void RunTests(); 
     
    4349        const Matrix<float> &GetEigenTransform() { return m_EigenTransform; } 
    4450        const FeatureVec &GetFeatures() { return m_Features; } 
     51        const Vector<float> &GetMean() { return m_Mean; } 
     52         
     53         
     54        void Load(FILE *f); 
     55        void Save(FILE *f); 
    4556         
    4657private:         
  • foam/trunk/vision/src/Vector.h

    r86 r90  
    2323 
    2424template<class T> 
     25bool feq(T a, T b, T t=0.001) 
     26{ 
     27        return fabs(a-b)<t; 
     28} 
     29 
     30template<class T> 
    2531class Vector 
    2632{ 
     
    5157        void Zero() { SetAll(0); } 
    5258        bool IsInf(); 
     59        T Mean(); 
    5360 
    5461        Vector &operator=(const Vector &other); 
    5562        Vector operator+(const Vector &other) const; 
    5663        Vector operator-(const Vector &other) const; 
     64        Vector operator+(T v) const; 
     65        Vector operator-(T v) const; 
    5766        Vector operator*(T v) const; 
    5867        Vector operator/(T v) const; 
    5968        Vector &operator+=(const Vector &other); 
    6069        Vector &operator-=(const Vector &other); 
     70        Vector &operator+=(T v); 
     71        Vector &operator-=(T v); 
    6172        Vector &operator*=(T v); 
    6273        Vector &operator/=(T v); 
     74         
     75        void Save(FILE *f); 
     76        void Load(FILE *f); 
    6377         
    6478        static void RunTests(); 
     
    170184 
    171185template<class T> 
     186Vector<T> Vector<T>::operator+(T v) const 
     187{        
     188        Vector<T> ret(m_Size); 
     189        for (unsigned int i=0; i<m_Size; i++) 
     190        { 
     191                ret[i]=(*this)[i]+v; 
     192        } 
     193        return ret; 
     194} 
     195 
     196template<class T> 
     197Vector<T> Vector<T>::operator-(T v) const 
     198{ 
     199        Vector<T> ret(m_Size); 
     200        for (unsigned int i=0; i<m_Size; i++) 
     201        { 
     202                ret[i]=(*this)[i]-v; 
     203        } 
     204        return ret; 
     205} 
     206 
     207template<class T> 
    172208Vector<T> Vector<T>::operator*(T v) const 
    173209{        
     
    194230Vector<T> &Vector<T>::operator+=(const Vector &other) 
    195231{ 
     232        assert(m_Size==other.m_Size); 
    196233        for (unsigned int i=0; i<m_Size; i++) 
    197234        { 
     
    204241Vector<T> &Vector<T>::operator-=(const Vector &other) 
    205242{ 
     243        assert(m_Size==other.m_Size); 
    206244        for (unsigned int i=0; i<m_Size; i++) 
    207245        { 
     
    212250 
    213251template<class T> 
     252Vector<T> &Vector<T>::operator+=(T v) 
     253{ 
     254        for (unsigned int i=0; i<m_Size; i++) 
     255        { 
     256                (*this)[i]+=v; 
     257        } 
     258        return *this; 
     259} 
     260 
     261template<class T> 
     262Vector<T> &Vector<T>::operator-=(T v) 
     263{ 
     264        for (unsigned int i=0; i<m_Size; i++) 
     265        { 
     266                (*this)[i]-=v; 
     267        } 
     268        return *this; 
     269} 
     270 
     271template<class T> 
    214272Vector<T> &Vector<T>::operator*=(T v) 
    215273{ 
     
    231289} 
    232290 
     291template<class T> 
     292T Vector<T>::Mean() 
     293{ 
     294        T acc=0; 
     295        for (unsigned int i=0; i<m_Size; i++) 
     296        { 
     297                acc+=(*this)[i]; 
     298        } 
     299        return acc/(float)m_Size; 
     300} 
     301 
     302template<class T> 
     303void Vector<T>::Save(FILE* f) 
     304{ 
     305        int version = 1;         
     306        fwrite(&version,sizeof(version),1,f); 
     307        fwrite(&m_Size,sizeof(m_Size),1,f); 
     308        fwrite(m_Data,sizeof(T)*m_Size,1,f); 
     309} 
     310 
     311template<class T> 
     312void Vector<T>::Load(FILE* f) 
     313{ 
     314        int version;     
     315        fread(&version,sizeof(version),1,f); 
     316        fread(&m_Size,sizeof(m_Size),1,f); 
     317        m_Data=new T[m_Size]; 
     318        fread(m_Data,sizeof(T)*m_Size,1,f);      
     319} 
    233320 
    234321template<class T> 
  • foam/trunk/vision/src/main.cpp

    r89 r90  
    3333double scale = 1; 
    3434 
    35 Image ti("data/spot-1.png"); 
    36 PCA pca(ti.RGB2GRAY().Scale(30,40).NumElements()); 
    37  
    38 void TestPCA() 
    39 { 
    40         ti=ti.RGB2GRAY().Scale(30,40); 
     35//int w=50; 
     36//int h=80; 
     37int w=20; 
     38int h=30; 
     39 
     40PCA pca(w*h); 
     41Vector<float> params(50); 
     42Image src("data/dave.png"); 
     43 
     44void Recalc() 
     45{ 
    4146        glob_t g; 
    42         glob("data/pca/*.png",GLOB_PERIOD,NULL,&g); 
     47         
     48        glob("data/spacek-large/*.png",GLOB_PERIOD,NULL,&g); 
    4349        for (unsigned int n=0; n<g.gl_pathc; n++) 
    4450        { 
     
    4753                Image im(path); 
    4854                //im.SubMean(); 
    49                 pca.AddFeature(im.RGB2GRAY().Scale(30,40).ToFloatVector()); 
     55                Vector<float> v(im.Scale(w,h).RGB2GRAY().ToFloatVector()); 
     56                v-=v.Mean(); 
     57                pca.AddFeature(v); 
    5058        } 
    5159        globfree (&g); 
    52         pca.Calculate(); 
     60         
     61        pca.Calculate();  
     62} 
     63 
     64void TestPCA() 
     65{ 
     66        //Recalc(); 
     67        //FILE *f=fopen("spacek-20x30.pca", "wb"); 
     68        //pca.Save(f); 
     69         
     70        FILE *f=fopen("spacek-20x30.pca", "rb"); 
     71        pca.Load(f); 
     72         
     73        fclose(f); 
     74         
     75        pca.Compress(0,50); 
     76        src = src.Scale(w,h).RGB2GRAY(); 
     77        Vector<float> d(src.ToFloatVector());    
     78        params=pca.Project(d);   
     79        params[0]=1; 
     80 
    5381} 
    5482 
     
    205233        ////////////////////////////////// 
    206234        // Matrix tests 
    207         Matrix<float>::RunTests(); 
     235        //Matrix<float>::RunTests(); 
    208236 
    209237        ////////////////////////////////// 
     
    263291        /////////////////////////////////// 
    264292        // PCA display 
    265          
    266293        camera.Clear(); 
    267294 
    268         for (unsigned int i=0; i<pca.GetFeatures().size(); i++) 
    269         { 
    270                 camera.Blit(Image(30,40,1,pca.GetFeatures()[i]),(i%20)*32,(i/20)*42); 
    271         } 
    272          
    273         for (unsigned int i=0; i<pca.GetEigenTransform().GetCols(); i++) 
    274         { 
    275                 camera.Blit(Image(30,40,1,pca.GetEigenTransform().GetColVector(i),5),(i%20)*32,200+(i/20)*42); 
    276         } 
    277          
    278         pca.GetEigenValues().Print(); 
     295        //for (unsigned int i=0; i<pca.GetFeatures().size(); i++) 
     296        //{ 
     297        //      camera.Blit(Image(30,40,1,pca.GetFeatures()[i]),(i%20)*32,(i/20)*42); 
     298        //} 
     299         
     300        for (unsigned int i=0; i<pca.GetEigenTransform().GetRows(); i++) 
     301        { 
     302                camera.Blit(Image(w,h,1,pca.GetEigenTransform().GetRowVector(i)*2+pca.GetMean()),(i%20)*(w+2),0+(i/20)*(h+2)); 
     303        } 
     304         
     305        camera.Blit(Image(w,h,1,pca.Synth(params)),200,200); 
     306        camera.Blit(src,100,200); 
    279307 
    280308    cvShowImage("result", camera.m_Image); 
Note: See TracChangeset for help on using the changeset viewer.