EMAN2
Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Attributes | Private Attributes
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:
maskmask image
nvecnumber of desired basis vectors
nimgtotal number of input images, required even with insert_image()

Definition at line 305 of file analyzer.h.


Constructor & Destructor Documentation

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

Definition at line 308 of file analyzer.h.

Referenced by NEW().

: 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 840 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, EMAN::EMData::set_attr(), EMAN::EMData::set_size(), V, and X.

{
// Allocate the working space
gsl_vector *work=gsl_vector_alloc(nimg);
gsl_vector *S=gsl_vector_alloc(nimg);
gsl_matrix *V=gsl_matrix_alloc(nimg,nimg);
gsl_matrix *X=gsl_matrix_alloc(nimg,nimg);


// Do the decomposition. All the real work is here
gsl_linalg_SV_decomp_mod (A,X, V, S, work);
//else gsl_linalg_SV_decomp_jacobi(A,V,S);

vector<EMData*> ret;
//unpack the results and write the output file
float *md=mask->get_data();
size_t totpix=mask->get_xsize()*mask->get_ysize()*mask->get_zsize();
for (int k=0; k<nvec; k++) {
        EMData *img = new EMData;
        img->set_size(mask->get_xsize(),mask->get_ysize(),mask->get_zsize());

        float  *d=img->get_data();
        for (size_t i=0,j=0; i<totpix; ++i) {
                if (md[i]) {
                        d[i]=(float)gsl_matrix_get(A,j,k);
                        j++;
                }
        }
        img->set_attr( "eigval", gsl_vector_get(S,k));
        ret.push_back(img);
}

gsl_vector_free(work);
gsl_vector_free(S);
gsl_matrix_free(V);
gsl_matrix_free(X);

gsl_matrix_free(A);
A=NULL;
mask=NULL;

return ret;
}
string EMAN::SVDAnalyzer::get_desc ( ) const [inline, virtual]

Get the Analyzer's description.

Returns:
The Analyzer's description.

Implements EMAN::Analyzer.

Definition at line 327 of file analyzer.h.

                {
                        return "Singular Value Decomposition from GSL. Comparable to pca";
                }
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 322 of file analyzer.h.

References NAME.

                {
                        return NAME;
                }
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 339 of file analyzer.h.

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

                {
                        TypeDict d;
                        d.put("mask", EMObject::EMDATA, "mask image");
                        d.put("nvec", EMObject::INT, "number of desired basis vectors");
                        d.put("nimg", EMObject::INT, "total number of input images, required even with insert_image()");
                        return d;
                }
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 818 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.

{
        if (mask==0)
                throw NullPointerException("Null mask image pointer, set_params() first");

        // count pixels under mask
        size_t totpix=mask->get_xsize()*mask->get_ysize()*mask->get_zsize();
        float  *d=image->get_data();
        float *md=mask ->get_data();
        for (size_t i=0,j=0; i<totpix; ++i) {
                if (md[i]) {
                        gsl_matrix_set(A,j,nsofar,d[i]);
                        j++;
                }
        }
        nsofar++;

   return 0;
}
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 312 of file analyzer.h.

References EMAN::Analyzer::images.

                                                                            {
                        vector<EMData*>::const_iterator iter;
                        for(iter=image_list.begin(); iter!=image_list.end(); ++iter) {
                                images.push_back(*iter);
                        }
                        return 0;
                }
static Analyzer* EMAN::SVDAnalyzer::NEW ( ) [inline, static]

Definition at line 332 of file analyzer.h.

References SVDAnalyzer().

                {
                        return new SVDAnalyzer();
                }
void SVDAnalyzer::set_params ( const Dict new_params) [virtual]

Set the Analyzer parameters using a key/value dictionary.

Parameters:
new_paramsA dictionary containing the new parameters.

Reimplemented from EMAN::Analyzer.

Definition at line 884 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.

{
        params = new_params;
        mask = params["mask"];
        nvec = params["nvec"];
        nimg = params["nimg"];

        // count pixels under mask
        pixels=0;
        size_t totpix=mask->get_xsize()*mask->get_ysize()*mask->get_zsize();
        float *d=mask->get_data();
        for (size_t i=0; i<totpix; ++i) if (d[i]) ++pixels;

        printf("%d,%d\n",pixels,nimg);
        A=gsl_matrix_alloc(pixels,nimg);
        nsofar=0;
}

Member Data Documentation

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

Definition at line 358 of file analyzer.h.

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

Definition at line 351 of file analyzer.h.

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

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

Definition at line 348 of file analyzer.h.

Referenced by get_name().

int EMAN::SVDAnalyzer::nimg [protected]

Definition at line 354 of file analyzer.h.

Referenced by analyze(), and set_params().

Definition at line 357 of file analyzer.h.

Referenced by insert_image(), and set_params().

int EMAN::SVDAnalyzer::nvec [protected]

Definition at line 352 of file analyzer.h.

Referenced by analyze(), and set_params().

int EMAN::SVDAnalyzer::pixels [protected]

Definition at line 353 of file analyzer.h.

Referenced by set_params().


The documentation for this class was generated from the following files: