1 // This is brl/bbas/brad/brad_eigenspace.h
2 #ifndef brad_eigenspace_h
3 #define brad_eigenspace_h
4 //:
5 // \file
6 // \brief Process image regions to produce an eigensystem from feature vectors
7 // \author J.L. Mundy
8 // \date June 30, 2011
9 //
10 // \verbatim
11 //  Modifications
12 //   <none yet>
13 // \endverbatim
14 //
15 //
16 //
17 // Each input resource is decomposed into a grid of blocks. A feature
18 // vector is extracted from the block, such as a histogram. A covariance
19 // matrix of the feature vectors is computed as well as the associated
20 // eigensystem. This class can be stored as binary to retain the
21 // eigenvector projection needed for classifying new images. The current
22 // implementation assumes that only the top three eigenvectors are needed
23 // by classifiers. This choice is convenient for display as color images
24 // and 3-d histograms.
25 //
26 // * nib - number of columns in a block
27 // * njb - number of rows in a block
28 // * funct - a functor class that computes the feature vector
29 //           the functor must return a vnl_vector<float> from
30 //           its (vil_image_view<float> const&) operator
31 //           and has a size() method to indicate the number of vector elements
32 //
33 #include <vector>
34 #include <iostream>
35 #include <vsl/vsl_binary_io.h>
36 #ifdef _MSC_VER
37 #  include <vcl_msvc_warnings.h>
38 #endif
39 #include <vnl/vnl_vector.h>
40 #include <vnl/vnl_matrix.h>
41 #include <vil/vil_image_resource.h>
42 #include <vil/vil_image_view.h>
43 #include <bsta/bsta_joint_histogram_3d.h>
44 #include <brad/brad_eigenspace_base.h>
45 
46 template <class T>
47 class brad_eigenspace : public brad_eigenspace_base
48 {
49  public:
brad_eigenspace()50   brad_eigenspace()
51   : nib_(0), njb_(0), covar_valid_(false), eigensystem_valid_(false) {}
52 
53   //: constructor with block size and functor
brad_eigenspace(unsigned nib,unsigned njb,T funct)54   brad_eigenspace(unsigned nib, unsigned njb, T funct)
55   : funct_(funct), nib_(nib), njb_(njb), covar_valid_(false), eigensystem_valid_(false) {}
56 
57   ~brad_eigenspace() override = default;
58 
feature_vector_type()59   std::string feature_vector_type() override {return funct_.type();}
60 
61   //: compute the covariance matrix from image resources
62   bool compute_covariance_matrix(std::vector<vil_image_resource_sptr> const& rescs);
63     //: compute the covariance matrix from a random fraction of image resources
64   bool compute_covariance_matrix_rand(std::vector<vil_image_resource_sptr> const& rescs, double frac = 1.0, unsigned nit=256, unsigned njt=256);
65   //: compute the covariance matrix using blocked image cache
66   bool compute_covariance_matrix_blocked(std::vector<vil_image_resource_sptr> const& rescs, unsigned nit=256, unsigned njt=256);
67   //: compute the eigensystem from the covariance matrix
68   bool compute_eigensystem();
69 
70   //: projection of image blocks onto the top three eigenvectors
71   bool compute_eigenimage(vil_image_resource_sptr const& resc,
72                           std::string const& output_path);
73 
74   //: projection of image blocks onto the top three eigenvectors, per pixel. assumes input and output can fit in memory
75   bool compute_eigenimage_pixel(vil_image_view<float> const& input,
76                                 vil_image_view<float>& eignimage);
77 
78   //: create and update an initial histogram of feature vectors projected onto the first three eigenvectors.Defines the min, max range of the histogram cells.
79   // note, does not update the histogram counts
80   bool init_histogram(vil_image_resource_sptr const& resc, unsigned nbins,
81                       bsta_joint_histogram_3d<float>& hist);
82 
83   //: update an existing histogram from an image resource
84   bool update_histogram(vil_image_resource_sptr const& resc,
85                         bsta_joint_histogram_3d<float>& hist);
86 
87   //: init from a set of images
88   bool init_histogram(std::vector<vil_image_resource_sptr> const& rescs,
89                       unsigned nbins,
90                       bsta_joint_histogram_3d<float>& hist);
91   //: update from a set of images
92   bool update_histogram(std::vector<vil_image_resource_sptr> const& rescs,
93                         bsta_joint_histogram_3d<float>& hist);
94 
95   //: init from a set of images, random selection of tiles
96   bool init_histogram_rand(std::vector<vil_image_resource_sptr> const& rescs,
97                            unsigned nbins,
98                            bsta_joint_histogram_3d<float>& hist,
99                            double frac = 1.0, unsigned nit=256,
100                            unsigned njt=256);
101 
102   //: update from a set of images, random selection of tiles
103   bool update_histogram_rand(std::vector<vil_image_resource_sptr> const& rescs,
104                              bsta_joint_histogram_3d<float>& hist,
105                              double frac = 1.0, unsigned nit=256,
106                              unsigned njt=256);
107     //: init from a set of images, using a blocked cache
108   bool init_histogram_blocked(std::vector<vil_image_resource_sptr> const& rescs,
109                               unsigned nbins,
110                               bsta_joint_histogram_3d<float>& hist,
111                               unsigned nit=256, unsigned njt=256);
112 
113 
114   //: update from a set of images, using a blocked cache
115   bool update_histogram_blocked(std::vector<vil_image_resource_sptr> const& rescs,
116                              bsta_joint_histogram_3d<float>& hist,
117                              unsigned nit=256, unsigned njt=256);
118   //: classify an image as to atmospheric content, nit, njt are block size
119   bool classify_image(vil_image_resource_sptr const& resc,
120                       bsta_joint_histogram_3d<float> const& no_atmos,
121                       bsta_joint_histogram_3d<float> const& atmos,
122                       unsigned nit, unsigned njt,
123                       std::string const& output_path);
124   bool classify_image(vil_image_resource_sptr const& resc,
125                       bsta_joint_histogram_3d<float> const& no_atmos,
126                       bsta_joint_histogram_3d<float> const& atmos,
127                       unsigned nit, unsigned njt,
128                       vil_image_resource_sptr& out_resc,
129                       vil_image_resource_sptr& out_resc_orig_size);
130 
131   //: classify an image as to atmospheric content, nit, njt are block size, assumes the resulting image can fit in memory. uses a sliding window to classify each pixel.
132   bool classify_image_pixel(vil_image_view<float> const& image,
133                             bsta_joint_histogram_3d<float> const& no_atmos,
134                             bsta_joint_histogram_3d<float> const& atmos,
135                             float prob_ratio,
136                             vil_image_view<float>& class_image);
137 
138   //: print sizes and other info (bool for calling consistency)
139   bool print(std::ostream& os = std::cout) const;
140 
141   //: accessors, setters
nib()142   unsigned nib() const {return nib_;}
njb()143   unsigned njb() const {return njb_;}
mean()144   vnl_vector<double> mean() const {return mean_;}
covariance()145   vnl_matrix<double> covariance() const {return covar_;}
eigenvalues()146   vnl_vector<double> eigenvalues() const {return eigenvalues_;}
eigenvectors()147   vnl_matrix<double> eigenvectors() const {return eigenvectors_;}
functor()148   T functor() const {return funct_;}
set_nib(unsigned nib)149   void set_nib(unsigned nib) {nib_=nib;}
set_njb(unsigned njb)150   void set_njb(unsigned njb) {njb_=njb;}
set_mean_covar(vnl_vector<double> const & mean,vnl_matrix<double> const & covar)151   void set_mean_covar(vnl_vector<double> const& mean,
152                       vnl_matrix<double> const& covar)
153   { mean_ = mean; covar_ = covar; covar_valid_ = true; }
set_eigensystem(vnl_vector<double> const & eigv,vnl_matrix<double> const & evecs)154   void set_eigensystem(vnl_vector<double> const& eigv,
155                        vnl_matrix<double> const& evecs)
156   { eigenvalues_=eigv; eigenvectors_=evecs; eigensystem_valid_ = true; }
set_functor(T funct)157   void set_functor(T funct) {funct_ = funct;}
158 
159  private:
160   //: An instance of the functor
161   T funct_;
162   //: number of cols in view block
163   unsigned nib_;
164   //: number of rows in view block
165   unsigned njb_;
166   //: feature vector mean
167   vnl_vector<double> mean_;
168   //: feature vector covariance matrix
169   vnl_matrix<double> covar_;
170   //: eigenvalues of the covariance matrix
171   vnl_vector<double> eigenvalues_;
172   //: eigenvectors of the covariance matrix
173   vnl_matrix<double> eigenvectors_;
174   //: flags
175   bool covar_valid_;
176   bool eigensystem_valid_;
177 };
178 
179 //: Binary save brad_eigenspace to stream.
180 template <class T>
181 void vsl_b_write(vsl_b_ostream &os, const brad_eigenspace<T>& ep);
182 
183 
184 //: Binary load brad_eigenspace from stream.
185 template <class T>
186 void vsl_b_read(vsl_b_istream &is, brad_eigenspace<T>& fv);
187 
188 
189 //: Print summary
190 template <class T>
191 void
192 vsl_print_summary(std::ostream &os, const brad_eigenspace<T>& fv);
193 
194 #include <brad/brad_eigenspace_sptr.h>
195 
196 #endif // brad_eigenspace_h
197