EMAN2
Public Member Functions | Static Public Member Functions | Static Public Attributes
EMAN::RT3DSphereAligner Class Reference

3D rotational and translational alignment using spherical sampling, can reduce the search space based on symmetry. More...

#include <aligner.h>

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

List of all members.

Public Member Functions

virtual EMDataalign (EMData *this_img, EMData *to_img, const string &cmp_name="sqeuclidean", const Dict &cmp_params=Dict()) const
 See Aligner comments for more details.
virtual EMDataalign (EMData *this_img, EMData *to_img) const
 See Aligner comments for more details.
virtual vector< Dictxform_align_nbest (EMData *this_img, EMData *to_img, const unsigned int nsoln, const string &cmp_name, const Dict &cmp_params) const
 See Aligner comments for more details.
virtual string get_name () const
 Get the Aligner's name.
virtual string get_desc () const
virtual TypeDict get_param_types () const

Static Public Member Functions

static AlignerNEW ()

Static Public Attributes

static const string NAME = "rotate_translate_3d"

Detailed Description

3D rotational and translational alignment using spherical sampling, can reduce the search space based on symmetry.

can also make use of different OrientationGenerators (random, for example) 2X more efficient than the RT3DGridAligner The aligner actually aligns the reference to the 'moving' and then takes the inverse of the resulting transform. This is necessary because, in the case of symmetry (i.e. not c1), the reference symmetry axis must be aligned to the EMAN2 symmetry axis, restricting the search space to the asymmetrical points on a sphere. We note that if the reference symmetry axis is not aligned to the EMAN2 symmetry axis, the best thing is to do a full search (i.e. specify sym='c1') unless you really know what you are doing!

Parameters:
symThe symmtery to use as the basis of the spherical sampling
orietgenAdvanced. The orientation generation strategy
deltaAngle the separates points on the sphere. This is exclusive of the 'n' paramater
nAn alternative to the delta argument, this is the number of points you want generated on the sphere
dphiThe angle increment in the phi direction
lphiLower bound for the phi direction
uphiUpper bound for the phi direction
dotransDo a translational search
searchThe maximum length of the detectable translational shift - if you supply this parameter you can not supply the maxshiftx, maxshifty or maxshiftz parameters. Each approach is mutually exclusive
searchxThe maximum length of the detectable translational shift in the x direction- if you supply this parameter you can not supply the maxshift parameters
searchyThe maximum length of the detectable translational shift in the y direction- if you supply this parameter you can not supply the maxshift parameters
searchzThe maximum length of the detectable translational shift in the z direction- if you supply this parameter you can not supply the maxshift parameters
verboseTurn this on to have useful information printed to standard out
Author:
John Flanagan and David Woolford
Date:
Feb 2011

Definition at line 1608 of file aligner.h.


Member Function Documentation

EMData * RT3DSphereAligner::align ( EMData this_img,
EMData to_img,
const string &  cmp_name = "sqeuclidean",
const Dict cmp_params = Dict() 
) const [virtual]

See Aligner comments for more details.

Implements EMAN::Aligner.

Definition at line 2612 of file aligner.cpp.

References EMAN::EMData::process(), EMAN::EMData::set_attr(), t, and EMAN::Aligner::xform_align_nbest().

Referenced by align().

{

        vector<Dict> alis = xform_align_nbest(this_img,to,1,cmp_name,cmp_params);

        Dict t;
        Transform* tr = (Transform*) alis[0]["xform.align3d"];
        t["transform"] = tr;
        EMData* soln = this_img->process("xform",t);
        soln->set_attr("xform.align3d",tr);
        delete tr; tr = 0;

        return soln;

}
virtual EMData* EMAN::RT3DSphereAligner::align ( EMData this_img,
EMData to_img 
) const [inline, virtual]

See Aligner comments for more details.

Implements EMAN::Aligner.

Definition at line 1617 of file aligner.h.

References align().

                        {
                                return align(this_img, to_img, "sqeuclidean", Dict());
                        }
virtual string EMAN::RT3DSphereAligner::get_desc ( ) const [inline, virtual]

Implements EMAN::Aligner.

Definition at line 1632 of file aligner.h.

                        {
                                return "3D rotational and translational alignment using spherical sampling. Can reduce the search space if symmetry is supplied";
                        }
virtual string EMAN::RT3DSphereAligner::get_name ( ) const [inline, virtual]

Get the Aligner's name.

Each Aligner is identified by a unique name.

Returns:
The Aligner's name.

Implements EMAN::Aligner.

Definition at line 1627 of file aligner.h.

References NAME.

                        {
                                return NAME;
                        }
virtual TypeDict EMAN::RT3DSphereAligner::get_param_types ( ) const [inline, virtual]

Implements EMAN::Aligner.

Definition at line 1642 of file aligner.h.

References EMAN::EMObject::BOOL, EMAN::EMObject::FLOAT, EMAN::EMObject::INT, EMAN::TypeDict::put(), EMAN::EMObject::STRING, and EMAN::EMObject::TRANSFORM.

                        {
                                TypeDict d;
                                d.put("sym", EMObject::STRING,"The symmtery to use as the basis of the spherical sampling. Default is c1 (asymmetry).");
                                d.put("orientgen", EMObject::STRING,"Advanced. The orientation generation strategy. Default is eman");
                                d.put("delta", EMObject::FLOAT,"Angle the separates points on the sphere. This is exclusive of the \'n\' paramater. Default is 10");
                                d.put("n", EMObject::INT,"An alternative to the delta argument, this is the number of points you want generated on the sphere. Default is OFF");
                                d.put("dphi", EMObject::FLOAT,"The angle increment in the phi direction. Default is 10");
                                d.put("phi0", EMObject::FLOAT,"Lower bound for phi. Default it 0");
                                d.put("phi1", EMObject::FLOAT,"Upper bound for phi. Default it 360");
                                d.put("dotrans", EMObject::BOOL,"Do a translational search. Default is True(1)");
                                d.put("search", EMObject::INT,"The maximum length of the detectable translational shift - if you supply this parameter you can not supply the maxshiftx, maxshifty or maxshiftz parameters. Each approach is mutually exclusive.");
                                d.put("searchx", EMObject::INT,"The maximum length of the detectable translational shift in the x direction- if you supply this parameter you can not supply the maxshift parameters. Default is 3.");
                                d.put("searchy", EMObject::INT,"The maximum length of the detectable translational shift in the y direction- if you supply this parameter you can not supply the maxshift parameters. Default is 3.");
                                d.put("searchz", EMObject::INT,"The maximum length of the detectable translational shift in the z direction- if you supply this parameter you can not supply the maxshift parameters. Default is 3");
                                d.put("initxform", EMObject::TRANSFORM,"The Transform storing the starting position. If unspecified the identity matrix is used");
                                d.put("verbose", EMObject::BOOL,"Turn this on to have useful information printed to standard out.");
                                return d;
                        }
static Aligner* EMAN::RT3DSphereAligner::NEW ( ) [inline, static]

Definition at line 1637 of file aligner.h.

                        {
                                return new RT3DSphereAligner();
                        }
vector< Dict > RT3DSphereAligner::xform_align_nbest ( EMData this_img,
EMData to_img,
const unsigned int  nsoln,
const string &  cmp_name,
const Dict cmp_params 
) const [virtual]

See Aligner comments for more details.

Reimplemented from EMAN::Aligner.

Definition at line 2628 of file aligner.cpp.

References EMAN::EMData::calc_ccf(), EMAN::EMData::calc_max_location_wrap(), calc_max_location_wrap_cuda(), EMAN::Cmp::cmp(), copy(), data, EMAN::EMData::do_fft(), EMAN::Symmetry3D::gen_orientations(), EMAN::EMData::get_ndim(), get_stats_cuda(), EMAN::EMData::get_value_at_wrap(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), EMAN::EMData::get_zsize(), EMAN::Dict::has_key(), ImageDimensionException, InvalidParameterException, EMAN::Transform::invert(), EMAN::Aligner::params, CudaPeakInfo::peak, phi, EMAN::EMData::process(), EMAN::EMData::process_inplace(), CudaPeakInfo::px, CudaPeakInfo::py, CudaPeakInfo::pz, EMAN::Dict::set_default(), EMAN::Transform::set_trans(), sqrt(), and t.

                                                                                                                                                                 {

        if ( this_img->get_ndim() != 3 || to->get_ndim() != 3 ) {
                throw ImageDimensionException("This aligner only works for 3D images");
        }

        int searchx = 0;
        int searchy = 0;
        int searchz = 0;
         
        bool dotrans = params.set_default("dotrans",1);
        if (params.has_key("search")) {
                vector<string> check;
                check.push_back("searchx");
                check.push_back("searchy");
                check.push_back("searchz");
                for(vector<string>::const_iterator cit = check.begin(); cit != check.end(); ++cit) {
                        if (params.has_key(*cit)) throw InvalidParameterException("The search parameter is mutually exclusive of the searchx, searchy, and searchz parameters");
                }
                int search  = params["search"];
                searchx = search;
                searchy = search;
                searchz = search;
        } else {
                searchx = params.set_default("searchx",3);
                searchy = params.set_default("searchy",3);
                searchz = params.set_default("searchz",3);
        }

        Transform* initxform;
        if (params.has_key("initxform") ) {
                // Unlike the 2d refine aligner, this class doesn't require the starting transform's
                // parameters to form the starting guess. Instead the Transform itself
                // is perturbed carefully (using quaternion rotation) to overcome problems that arise
                // when you use orthogonally-based Euler angles
                initxform = params["initxform"];
        }else {
                initxform = new Transform(); // is the identity
        }
        
        float lphi = params.set_default("phi0",0.0f);
        float uphi = params.set_default("phi1",360.0f);
        float dphi = params.set_default("dphi",10.f);
        float threshold = params.set_default("threshold",0.f);
        if (threshold < 0.0f) throw InvalidParameterException("The threshold parameter must be greater than or equal to zero");
        bool verbose = params.set_default("verbose",false);
        
        //in case we arre aligning tomos
        Dict altered_cmp_params(cmp_params);
        if (cmp_name == "ccc.tomo") {
                altered_cmp_params.set_default("searchx", searchx);
                altered_cmp_params.set_default("searchy", searchy);
                altered_cmp_params.set_default("searchz", searchz);
                altered_cmp_params.set_default("norm", true);
        }

        vector<Dict> solns;
        if (nsoln == 0) return solns; // What was the user thinking?
        for (unsigned int i = 0; i < nsoln; ++i ) {
                Dict d;
                d["score"] = 1.e24;
                Transform t; // identity by default
                d["xform.align3d"] = &t; // deep copy is going on here
                solns.push_back(d);
        }

        Dict d;
        d["inc_mirror"] = true; // This should probably always be true for 3D case. If it ever changes we might have to make inc_mirror a parameter
        if ( params.has_key("delta") && params.has_key("n") ) {
                throw InvalidParameterException("The delta and n parameters are mutually exclusive in the RT3DSphereAligner aligner");
        } else if ( params.has_key("n") ) {
                d["n"] = params["n"];
        } else {
                // If they didn't specify either then grab the default delta - if they did supply delta we're still safe doing this
                d["delta"] = params.set_default("delta",10.f);
        }

        if ((string)params.set_default("orientgen","eman")=="eman") d["perturb"]=0;
        Symmetry3D* sym = Factory<Symmetry3D>::get((string)params.set_default("sym","c1"));
        vector<Transform> transforms = sym->gen_orientations((string)params.set_default("orientgen","eman"),d);

        bool tomography = (cmp_name == "ccc.tomo") ? 1 : 0;
        
        //precompute fixed FT, saves a LOT of time!!!
        EMData * this_imgfft = 0;
        if(dotrans || tomography){
                this_imgfft = this_img->do_fft();
        }
        
#ifdef EMAN2_USING_CUDA 
        if(EMData::usecuda == 1) {
                cout << "Using CUDA for 3D alignment" << endl;
                if(!to->getcudarodata()) to->copy_to_cudaro(); // Safer call
                if(!this_img->getcudarwdata()) this_img->copy_to_cuda();
                if(this_imgfft) this_imgfft->copy_to_cuda();
        }
#endif

        Transform trans = Transform();
        Cmp* c = Factory <Cmp>::get(cmp_name, cmp_params);
        
        bool use_cpu = true;
        for(vector<Transform>::const_iterator trans_it = transforms.begin(); trans_it != transforms.end(); trans_it++) {
                Dict params = trans_it->get_params("eman");
                
                if (verbose) {
                        float alt = params["alt"];
                        float az = params["az"];
                        cout << "Trying angle alt: " << alt << " az: " << az << endl;
                }

                for( float phi = lphi; phi < uphi; phi += dphi ) { 
                        params["phi"] = phi;
                        Transform t(params);
                        t = t*(*initxform);
                        
                        EMData* transformed;
                        transformed = to->process("xform",Dict("transform",&t));
                                
                        //need to do things a bit diffrent if we want to compare two tomos
                        float best_score = 0.0f;
                        // Dotrans is effectievly ignored for tomography
                        if(dotrans || tomography){
                                EMData* ccf = transformed->calc_ccf(this_imgfft);
#ifdef EMAN2_USING_CUDA 
                                if(EMData::usecuda == 1){
                                        // I use the following code rather than ccc.tomo to avoid doing two CCCs
                                        use_cpu = false;
                                        CudaPeakInfo* data = calc_max_location_wrap_cuda(ccf->getcudarwdata(), ccf->get_xsize(), ccf->get_ysize(), ccf->get_zsize(), searchx, searchy, searchz);
                                        trans.set_trans((float)-data->px, (float)-data->py, (float)-data->pz);
                                        t = trans*t;    //composite transform to reflect the fact that we have done a rotation first and THEN a transformation
                                        if (tomography) {
                                                float2 stats = get_stats_cuda(ccf->getcudarwdata(), ccf->get_xsize(), ccf->get_ysize(), ccf->get_zsize());
                                                best_score = -(data->peak - stats.x)/sqrt(stats.y); // Normalize, this is better than calling the norm processor since we only need to normalize one point
                                        } else {
                                                best_score = -data->peak;
                                        }
                                        delete data;
                                }
#endif
                                if(use_cpu){
                                        // I use the following code rather than ccc.tomo to avoid doing two CCCs
                                        if(tomography) ccf->process_inplace("normalize");
                                        IntPoint point = ccf->calc_max_location_wrap(searchx,searchy,searchz);
                                        trans.set_trans((float)-point[0], (float)-point[1], (float)-point[2]);
                                        t = trans*t;    //composite transform to reflect the fact that we have done a rotation first and THEN a transformation
                                        best_score = -ccf->get_value_at_wrap(point[0], point[1], point[2]);
                                }
                                delete ccf; ccf =0;
                                delete transformed; transformed = 0;// this is to stop a mem leak
                        }

                        if(!tomography){
                                if(!transformed) transformed = to->process("xform",Dict("transform",&t));
                                best_score = c->cmp(this_img,transformed);
                                delete transformed; transformed = 0;
                        }

                        unsigned int j = 0;
                        //cout << "alt " <<float(t.get_rotation("eman").get("alt")) << " az " << float(t.get_rotation("eman").get("az")) << " phi " << float(t.get_rotation("eman").get("phi")) << endl;
                        for ( vector<Dict>::iterator it = solns.begin(); it != solns.end(); ++it, ++j ) {
                                if ( (float)(*it)["score"] > best_score ) { // Note greater than - EMAN2 preferes minimums as a matter of policy
                                        vector<Dict>::reverse_iterator rit = solns.rbegin();
                                        copy(rit+1,solns.rend()-j,rit);
                                        Dict& d = (*it);
                                        d["score"] = best_score;
                                        t.invert(); //We actually moved the ref onto the moving, so we need to invert to do the opposite(this is done b/c the ref is aligned to the sym axis, whereas the mvoing is not)
                                        d["xform.align3d"] = &t; // deep copy is going on here
                                        break;
                                }
                        }

                }
        }
        
        if(this_imgfft) {delete this_imgfft; this_imgfft = 0;}
        if(sym!=0) delete sym;
        if (c != 0) delete c;
        
        return solns;

}

Member Data Documentation

const string RT3DSphereAligner::NAME = "rotate_translate_3d" [static]

Definition at line 1662 of file aligner.h.

Referenced by get_name().


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