analyzer.h

Go to the documentation of this file.
00001 
00005 /*
00006  * Author: Steven Ludtke, 04/10/2003 (sludtke@bcm.edu)
00007  * Copyright (c) 2000-2006 Baylor College of Medicine
00008  *
00009  * This software is issued under a joint BSD/GNU license. You may use the
00010  * source code in this file under either license. However, note that the
00011  * complete EMAN2 and SPARX software packages have some GPL dependencies,
00012  * so you are responsible for compliance with the licenses of these packages
00013  * if you opt to use BSD licensing. The warranty disclaimer below holds
00014  * in either instance.
00015  *
00016  * This complete copyright notice must be included in any revised version of the
00017  * source code. Additional authorship citations may be added, but existing
00018  * author citations must be preserved.
00019  *
00020  * This program is free software; you can redistribute it and/or modify
00021  * it under the terms of the GNU General Public License as published by
00022  * the Free Software Foundation; either version 2 of the License, or
00023  * (at your option) any later version.
00024  *
00025  * This program is distributed in the hope that it will be useful,
00026  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00027  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
00028  * GNU General Public License for more details.
00029  *
00030  * You should have received a copy of the GNU General Public License
00031  * along with this program; if not, write to the Free Software
00032  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
00033  *
00034  * */
00035 
00036 #ifndef eman_analyzer_h__
00037 #define eman_analyzer_h__
00038 
00039 #include "emobject.h"
00040 #include <gsl/gsl_linalg.h>
00041 using std::vector;
00042 
00043 namespace EMAN
00044 {
00045         class EMData;
00046 
00060         class Analyzer
00061         {
00062           public:
00063                 Analyzer() {}
00064 
00065                 virtual ~Analyzer()
00066                 {}
00067 
00072                 virtual int insert_image(EMData * image) = 0;
00073 
00078                 virtual int insert_images_list(vector<EMData *> image_list);
00079 
00083                 virtual vector<EMData*> analyze() = 0;
00084 
00088                 virtual string get_name() const = 0;
00089 
00093                 virtual string get_desc() const = 0;
00094 
00098                 virtual void set_params(const Dict & new_params)
00099                 {
00100                         params = new_params;
00101                 }
00102 
00106                 virtual Dict get_params() const
00107                 {
00108                         return params;
00109                 }
00110 
00117                 virtual TypeDict get_param_types() const = 0;
00118 
00119           protected:
00120                 mutable Dict params;
00121                 vector<EMData *> images;
00122         };
00123 
00138         class KMeansAnalyzer:public Analyzer
00139         {
00140           public:
00141                 KMeansAnalyzer() : ncls(0),verbose(0),minchange(0),maxiter(100),mininclass(2),slowseed(0) {}
00142 
00143                 virtual int insert_image(EMData *image) {
00144                         images.push_back(image);
00145                         return 0;
00146                 }
00147 
00148                 virtual vector<EMData*> analyze();
00149 
00150                 string get_name() const
00151                 {
00152                         return NAME;
00153                 }
00154 
00155                 string get_desc() const
00156                 {
00157                         return "k-means classification";
00158                 }
00159 
00160                 static Analyzer * NEW()
00161                 {
00162                         return new KMeansAnalyzer();
00163                 }
00164 
00165                 void set_params(const Dict & new_params);
00166 
00167                 TypeDict get_param_types() const
00168                 {
00169                         TypeDict d;
00170                         d.put("verbose", EMObject::INT, "Display progress if set, more detail with larger numbers (9 max)");
00171                         d.put("ncls", EMObject::INT, "number of desired classes");
00172                         d.put("maxiter", EMObject::INT, "maximum number of iterations");
00173                         d.put("minchange", EMObject::INT, "Terminate if fewer than minchange members move in an iteration");
00174                         d.put("mininclass", EMObject::INT, "Minumum number of particles to keep a class as good (not enforced at termination");
00175                         d.put("slowseed",EMObject::INT, "Instead of seeding all classes at once, it will gradually increase the number of classes by adding new seeds in groups with large standard deviations");
00176                         d.put("calcsigmamean",EMObject::INT, "Computes standard deviation of the mean image for each class-average (center), and returns them at the end of the list of centers");
00177                         return d;
00178                 }
00179 
00180                 static const string NAME;
00181                 
00182           protected:
00183                 void update_centers(int sigmas=0);
00184                 void reclassify();
00185                 void reseed();
00186 
00187                 vector<EMData *> centers;
00188                 int ncls;       //number of desired classes
00189                 int verbose;
00190                 int minchange;
00191                 int maxiter;
00192                 int mininclass;
00193                 int nchanged;
00194                 int slowseed;
00195                 int calcsigmamean;
00196 
00197         };
00198 
00204         class SVDAnalyzer : public Analyzer
00205         {
00206           public:
00207                 SVDAnalyzer() : mask(0), nvec(0), nimg(0), A(NULL) {}
00208 
00209                 virtual int insert_image(EMData * image);
00210 
00211                 virtual int insert_images_list(vector<EMData *> image_list) {
00212                         vector<EMData*>::const_iterator iter;
00213                         for(iter=image_list.begin(); iter!=image_list.end(); ++iter) {
00214                                 images.push_back(*iter);
00215                         }
00216                         return 0;
00217                 }
00218 
00219                 virtual vector<EMData*> analyze();
00220 
00221                 string get_name() const
00222                 {
00223                         return NAME;
00224                 }
00225 
00226                 string get_desc() const
00227                 {
00228                         return "Singular Value Decomposition from GSL. Comparable to pca";
00229                 }
00230 
00231                 static Analyzer * NEW()
00232                 {
00233                         return new SVDAnalyzer();
00234                 }
00235 
00236                 void set_params(const Dict & new_params);
00237 
00238                 TypeDict get_param_types() const
00239                 {
00240                         TypeDict d;
00241                         d.put("mask", EMObject::EMDATA, "mask image");
00242                         d.put("nvec", EMObject::INT, "number of desired basis vectors");
00243                         d.put("nimg", EMObject::INT, "total number of input images, required even with insert_image()");
00244                         return d;
00245                 }
00246 
00247                 static const string NAME;
00248                 
00249           protected:
00250                 EMData * mask;
00251                 int nvec;       //number of desired principal components
00252                 int pixels;     // pixels under the mask
00253                 int nimg; // number of input images
00254 
00255                 private:
00256                 int nsofar;
00257                 gsl_matrix *A;
00258         };
00259 
00260         template <> Factory < Analyzer >::Factory();
00261 
00262         void dump_analyzers();
00263         map<string, vector<string> > dump_analyzers_list();
00264 }
00265 
00266 #endif  //eman_analyzer_h__

Generated on Tue Jun 11 12:40:21 2013 for EMAN2 by  doxygen 1.4.7