EMAN::SVDAnalyzer Class Reference

Singular Value Decomposition from GSL. More...

#include <analyzer.h>

Inheritance diagram for EMAN::SVDAnalyzer:

Inheritance graph
[legend]
Collaboration diagram for EMAN::SVDAnalyzer:

Collaboration graph
[legend]
List of all members.

Public Member Functions

 SVDAnalyzer ()
virtual int insert_image (EMData *image)
 insert a image to the list of input images
virtual int insert_images_list (vector< EMData * > image_list)
 insert a list of images to the list of input images
virtual vector< EMData * > analyze ()
 main function for Analyzer, analyze input images and create output images
string get_name () const
 Get the Analyzer's name.
string get_desc () const
 Get the Analyzer's description.
void set_params (const Dict &new_params)
 Set the Analyzer parameters using a key/value dictionary.
TypeDict get_param_types () const
 Get Analyzer parameter information in a dictionary.

Static Public Member Functions

static AnalyzerNEW ()

Static Public Attributes

static const string NAME = "svd_gsl"

Protected Attributes

EMDatamask
int nvec
int pixels
int nimg

Private Attributes

int nsofar
gsl_matrix * A

Detailed Description

Singular Value Decomposition from GSL.

Comparable to pca

Parameters:
mask mask image
nvec number of desired basis vectors
nimg total number of input images, required even with insert_image()

Definition at line 204 of file analyzer.h.


Constructor & Destructor Documentation

EMAN::SVDAnalyzer::SVDAnalyzer (  )  [inline]

Definition at line 207 of file analyzer.h.

Referenced by NEW().

00207 : mask(0), nvec(0), nimg(0), A(NULL) {}


Member Function Documentation

vector< EMData * > SVDAnalyzer::analyze (  )  [virtual]

main function for Analyzer, analyze input images and create output images

Returns:
vector<EMData *> result os images analysis

Implements EMAN::Analyzer.

Definition at line 753 of file analyzer.cpp.

References A, EMAN::EMData::get_data(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), EMAN::EMData::get_zsize(), mask, nimg, nvec, V, and X.

00754 {
00755 // Allocate the working space
00756 gsl_vector *work=gsl_vector_alloc(nimg);
00757 gsl_vector *S=gsl_vector_alloc(nimg);
00758 gsl_matrix *V=gsl_matrix_alloc(nimg,nimg);
00759 gsl_matrix *X=gsl_matrix_alloc(nimg,nimg);
00760 
00761 
00762 // Do the decomposition. All the real work is here
00763 gsl_linalg_SV_decomp_mod (A,X, V, S, work);
00764 //else gsl_linalg_SV_decomp_jacobi(A,V,S);
00765 
00766 vector<EMData*> ret;
00767 //unpack the results and write the output file
00768 float *md=mask->get_data();
00769 size_t totpix=mask->get_xsize()*mask->get_ysize()*mask->get_zsize();
00770 for (int k=0; k<nvec; k++) {
00771         EMData *img = new EMData;
00772         img->set_size(mask->get_xsize(),mask->get_ysize(),mask->get_zsize());
00773 
00774         float  *d=img->get_data();
00775         for (size_t i=0,j=0; i<totpix; ++i) {
00776                 if (md[i]) {
00777                         d[i]=(float)gsl_matrix_get(A,j,k);
00778                         j++;
00779                 }
00780         }
00781         img->set_attr( "eigval", gsl_vector_get(S,k));
00782         ret.push_back(img);
00783 }
00784 
00785 gsl_vector_free(work);
00786 gsl_vector_free(S);
00787 gsl_matrix_free(V);
00788 gsl_matrix_free(X);
00789 
00790 gsl_matrix_free(A);
00791 A=NULL;
00792 mask=NULL;
00793 
00794 return ret;
00795 }

string EMAN::SVDAnalyzer::get_desc (  )  const [inline, virtual]

Get the Analyzer's description.

Returns:
The Analyzer's description.

Implements EMAN::Analyzer.

Definition at line 226 of file analyzer.h.

00227                 {
00228                         return "Singular Value Decomposition from GSL. Comparable to pca";
00229                 }

string EMAN::SVDAnalyzer::get_name (  )  const [inline, virtual]

Get the Analyzer's name.

Each Analyzer is identified by a unique name.

Returns:
The Analyzer's name.

Implements EMAN::Analyzer.

Definition at line 221 of file analyzer.h.

References NAME.

00222                 {
00223                         return NAME;
00224                 }

TypeDict EMAN::SVDAnalyzer::get_param_types (  )  const [inline, virtual]

Get Analyzer parameter information in a dictionary.

Each parameter has one record in the dictionary. Each record contains its name, data-type, and description.

Returns:
A dictionary containing the parameter info.

Implements EMAN::Analyzer.

Definition at line 238 of file analyzer.h.

References EMAN::EMObject::EMDATA, EMAN::EMObject::INT, and EMAN::TypeDict::put().

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                 }

int SVDAnalyzer::insert_image ( EMData image  )  [virtual]

insert a image to the list of input images

Parameters:
image 
Returns:
int 0 for success, <0 for fail

Implements EMAN::Analyzer.

Definition at line 731 of file analyzer.cpp.

References A, EMAN::EMData::get_data(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), EMAN::EMData::get_zsize(), mask, nsofar, and NullPointerException.

00732 {
00733         if (mask==0)
00734                 throw NullPointerException("Null mask image pointer, set_params() first");
00735 
00736         // count pixels under mask
00737         size_t totpix=mask->get_xsize()*mask->get_ysize()*mask->get_zsize();
00738         float  *d=image->get_data();
00739         float *md=mask ->get_data();
00740         for (size_t i=0,j=0; i<totpix; ++i) {
00741                 if (md[i]) {
00742                         gsl_matrix_set(A,j,nsofar,d[i]);
00743                         j++;
00744                 }
00745         }
00746         nsofar++;
00747 
00748    return 0;
00749 }

virtual int EMAN::SVDAnalyzer::insert_images_list ( vector< EMData * >  image_list  )  [inline, virtual]

insert a list of images to the list of input images

Parameters:
image_list 
Returns:
int 0 for success, <0 for fail

Reimplemented from EMAN::Analyzer.

Definition at line 211 of file analyzer.h.

References EMAN::Analyzer::images.

00211                                                                             {
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                 }

static Analyzer* EMAN::SVDAnalyzer::NEW (  )  [inline, static]

Definition at line 231 of file analyzer.h.

References SVDAnalyzer().

00232                 {
00233                         return new SVDAnalyzer();
00234                 }

void SVDAnalyzer::set_params ( const Dict new_params  )  [virtual]

Set the Analyzer parameters using a key/value dictionary.

Parameters:
new_params A dictionary containing the new parameters.

Reimplemented from EMAN::Analyzer.

Definition at line 797 of file analyzer.cpp.

References A, EMAN::EMData::get_data(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), EMAN::EMData::get_zsize(), mask, nimg, nsofar, nvec, EMAN::Analyzer::params, and pixels.

00798 {
00799         params = new_params;
00800         mask = params["mask"];
00801         nvec = params["nvec"];
00802         nimg = params["nimg"];
00803 
00804         // count pixels under mask
00805         pixels=0;
00806         size_t totpix=mask->get_xsize()*mask->get_ysize()*mask->get_zsize();
00807         float *d=mask->get_data();
00808         for (size_t i=0; i<totpix; ++i) if (d[i]) ++pixels;
00809 
00810         printf("%d,%d\n",pixels,nimg);
00811         A=gsl_matrix_alloc(pixels,nimg);
00812         nsofar=0;
00813 }


Member Data Documentation

gsl_matrix* EMAN::SVDAnalyzer::A [private]

Definition at line 257 of file analyzer.h.

Referenced by analyze(), insert_image(), and set_params().

EMData* EMAN::SVDAnalyzer::mask [protected]

Definition at line 250 of file analyzer.h.

Referenced by analyze(), insert_image(), and set_params().

const string EMAN::SVDAnalyzer::NAME = "svd_gsl" [static]

Definition at line 247 of file analyzer.h.

Referenced by get_name().

int EMAN::SVDAnalyzer::nimg [protected]

Definition at line 253 of file analyzer.h.

Referenced by analyze(), and set_params().

int EMAN::SVDAnalyzer::nsofar [private]

Definition at line 256 of file analyzer.h.

Referenced by insert_image(), and set_params().

int EMAN::SVDAnalyzer::nvec [protected]

Definition at line 251 of file analyzer.h.

Referenced by analyze(), and set_params().

int EMAN::SVDAnalyzer::pixels [protected]

Definition at line 252 of file analyzer.h.

Referenced by set_params().


The documentation for this class was generated from the following files:
Generated on Tue Jun 11 12:42:50 2013 for EMAN2 by  doxygen 1.4.7