EMAN2
Defines | Functions
aligner.cpp File Reference
#include "emfft.h"
#include "cmp.h"
#include "aligner.h"
#include "averager.h"
#include "emdata.h"
#include "processor.h"
#include "util.h"
#include "symmetry.h"
#include <gsl/gsl_multimin.h>
#include "plugins/aligner_template.h"
Include dependency graph for aligner.cpp:

Go to the source code of this file.

Defines

#define EMAN2_ALIGNER_DEBUG   0

Functions

static double refalifn (const gsl_vector *v, void *params)
static void refalidf (const gsl_vector *v, void *params, gsl_vector *df)
static void refalifdf (const gsl_vector *v, void *params, double *f, gsl_vector *df)
static double refalifnfast (const gsl_vector *v, void *params)
static Transform refalin3d_perturbquat (const Transform *const t, const float &spincoeff, const float &n0, const float &n1, const float &n2, const float &x, const float &y, const float &z)
static double symquat (const gsl_vector *v, void *params)
static double refalifn3dquat (const gsl_vector *v, void *params)

Define Documentation

#define EMAN2_ALIGNER_DEBUG   0
Id:
aligner.cpp,v 1.276 2014/02/19 20:26:42 vern Exp

Definition at line 55 of file aligner.cpp.


Function Documentation

static void refalidf ( const gsl_vector *  v,
void *  params,
gsl_vector *  df 
) [static]

Definition at line 1617 of file aligner.cpp.

References refalifn().

Referenced by EMAN::RefineAlignerCG::align().

                                                                         {
        // we do this using a simple local difference estimate due to the expense of the calculation. 
        // The step has to be large enough for the similarity metric
        // To provide an accurate change in value. 
        static double lstep[4] = { 0.05, 0.05, 0.1, 0.01 }; 
        
        gsl_vector *vc = gsl_vector_alloc(v->size);
        gsl_vector_memcpy(vc,v);
        
        double f = refalifn(v,params);
        for (unsigned int i=0; i<v->size; i++) {
                // double *vp = gsl_vector_ptr(vc,i);
                double vp = gsl_vector_get(vc,i);
                // *vp+=lstep[i];
                gsl_vector_set(vc,i,vp+lstep[i]);
                double f2 = refalifn(vc,params);
                // *vp-=lstep[i];
                gsl_vector_set(vc,i,vp);
                
                gsl_vector_set(df,i,(f2-f)/lstep[i]);
        }
        
        gsl_vector_free(vc);
        return;
}
static void refalifdf ( const gsl_vector *  v,
void *  params,
double *  f,
gsl_vector *  df 
) [static]

Definition at line 1643 of file aligner.cpp.

References refalifn().

Referenced by EMAN::RefineAlignerCG::align().

                                                                                       {
        // we do this using a simple local difference estimate due to the expense of the calculation. 
        // The step has to be large enough for the similarity metric
        // To provide an accurate change in value. 
        static double lstep[4] = { 0.05, 0.05, 0.1, 0.01 }; 
        
        gsl_vector *vc = gsl_vector_alloc(v->size);
        gsl_vector_memcpy(vc,v);
        
        *f = refalifn(v,params);
        for (unsigned int i=0; i<v->size; i++) {
                // double *vp = gsl_vector_ptr(vc,i);
                double vp = gsl_vector_get(vc,i);
                // *vp+=lstep[i];
                gsl_vector_set(vc,i,vp+lstep[i]);
                double f2 = refalifn(vc,params);
                // *vp-=lstep[i];
                gsl_vector_set(vc,i,vp);
                
                gsl_vector_set(df,i,(f2-*f)/lstep[i]);
        }
        
        gsl_vector_free(vc);
        return;
}
static double refalifn ( const gsl_vector *  v,
void *  params 
) [static]

Definition at line 1585 of file aligner.cpp.

References EMAN::Cmp::cmp(), EMAN::Dict::has_key(), EMAN::EMData::mult(), EMAN::EMData::process(), EMAN::Transform::set_mirror(), EMAN::Transform::set_scale(), EMAN::Transform::set_trans(), t, x, and y.

Referenced by EMAN::RefineAlignerCG::align(), EMAN::RefineAligner::align(), refalidf(), and refalifdf().

{
        Dict *dict = (Dict *) params;

        double x = gsl_vector_get(v, 0);
        double y = gsl_vector_get(v, 1);
        double a = gsl_vector_get(v, 2);

        EMData *this_img = (*dict)["this"];
        EMData *with = (*dict)["with"];
        bool mirror = (*dict)["mirror"];

        Transform t(Dict("type","2d","alpha",static_cast<float>(a)));
        t.set_trans((float)x,(float)y);
        t.set_mirror(mirror);
        if (v->size>3) {
                float sca=(float)gsl_vector_get(v, 3);
                if (sca<.7 || sca>1.3) return 1.0e20;
                t.set_scale((float)gsl_vector_get(v, 3));
        }
        EMData *tmp = this_img->process("xform",Dict("transform",&t));
        if (dict->has_key("mask")) tmp->mult(*(EMData *)((*dict)["mask"]));

//      printf("GSL %f %f %f %d %f\n",x,y,a,mirror,(float)gsl_vector_get(v, 3));
        Cmp* c = (Cmp*) ((void*)(*dict)["cmp"]);
        double result = c->cmp(tmp,with);

        if (tmp != 0) delete tmp;
        
        return result;
}
static double refalifn3dquat ( const gsl_vector *  v,
void *  params 
) [static]

Definition at line 2001 of file aligner.cpp.

References EMAN::Cmp::cmp(), EMAN::EMData::process(), refalin3d_perturbquat(), t, x, and y.

Referenced by EMAN::Refine3DAlignerQuaternion::align().

{
        Dict *dict = (Dict *) params;

        double n0 = gsl_vector_get(v, 0);
        double n1 = gsl_vector_get(v, 1);
        double n2 = gsl_vector_get(v, 2);
        double x = gsl_vector_get(v, 3);
        double y = gsl_vector_get(v, 4);
        double z = gsl_vector_get(v, 5);

        EMData *this_img = (*dict)["this"];
        EMData *with = (*dict)["with"];

        Transform* t = (*dict)["transform"];
        float spincoeff = (*dict)["spincoeff"];

        Transform soln = refalin3d_perturbquat(t,spincoeff,(float)n0,(float)n1,(float)n2,(float)x,(float)y,(float)z);

        EMData *tmp = this_img->process("xform",Dict("transform",&soln));
        Cmp* c = (Cmp*) ((void*)(*dict)["cmp"]);
        double result = c->cmp(tmp,with);
        if ( tmp != 0 ) delete tmp;

        //cout << result << endl;
        return result;
}
static double refalifnfast ( const gsl_vector *  v,
void *  params 
) [static]

Definition at line 1669 of file aligner.cpp.

References EMAN::EMData::dot_rotate_translate(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), x, and y.

Referenced by EMAN::RefineAlignerCG::align(), and EMAN::RefineAligner::align().

{
        Dict *dict = (Dict *) params;
        EMData *this_img = (*dict)["this"];
        EMData *img_to = (*dict)["with"];
        bool mirror = (*dict)["mirror"];

        double x = gsl_vector_get(v, 0);
        double y = gsl_vector_get(v, 1);
        double a = gsl_vector_get(v, 2);

        double r = this_img->dot_rotate_translate(img_to, (float)x, (float)y, (float)a, mirror);
        int nsec = this_img->get_xsize() * this_img->get_ysize();
        double result = 1.0 - r / nsec;

//      cout << result << " x " << x << " y " << y << " az " << a <<  endl;
        return result;
}
static Transform refalin3d_perturbquat ( const Transform *const  t,
const float &  spincoeff,
const float &  n0,
const float &  n1,
const float &  n2,
const float &  x,
const float &  y,
const float &  z 
) [static]

Definition at line 1951 of file aligner.cpp.

References EMAN::Vec3< Type >::normalize(), q, EMAN::Transform::set_trans(), and sqrt().

Referenced by EMAN::Refine3DAlignerQuaternion::align(), EMAN::SymAlignProcessorQuat::align(), refalifn3dquat(), and symquat().

{
        Vec3f normal(n0,n1,n2);
        normal.normalize();
        
        float omega = spincoeff*sqrt(n0*n0 + n1*n1 + n2*n2); // Here we compute the spin by the rotation axis vector length
        Dict d;
        d["type"] = "spin";
        d["omega"] = omega;
        d["n1"] = normal[0];
        d["n2"] = normal[1];
        d["n3"] = normal[2];
        //cout << omega << " " << normal[0] << " " << normal[1] << " " << normal[2] << " " << n0 << " " << n1 << " " << n2 << endl;
        
        Transform q(d);
        q.set_trans((float)x,(float)y,(float)z);
        
        q = q*(*t); //compose transforms        
        
        return q;
}
static double symquat ( const gsl_vector *  v,
void *  params 
) [static]

Definition at line 1973 of file aligner.cpp.

References EMAN::Cmp::cmp(), EMAN::EMData::process(), refalin3d_perturbquat(), t, x, and y.

Referenced by EMAN::SymAlignProcessorQuat::align().

{
        Dict *dict = (Dict *) params;

        double n0 = gsl_vector_get(v, 0);
        double n1 = gsl_vector_get(v, 1);
        double n2 = gsl_vector_get(v, 2);
        double x = gsl_vector_get(v, 3);
        double y = gsl_vector_get(v, 4);
        double z = gsl_vector_get(v, 5);

        EMData* volume = (*dict)["volume"];
        float spincoeff = (*dict)["spincoeff"];
        Transform* t = (*dict)["transform"];

        Transform soln = refalin3d_perturbquat(t,spincoeff,(float)n0,(float)n1,(float)n2,(float)x,(float)y,(float)z);

        EMData *tmp = volume->process("xform",Dict("transform",&soln));
        EMData *symtmp = tmp->process("xform.applysym",Dict("sym",(*dict)["sym"]));
        Cmp* c = (Cmp*) ((void*)(*dict)["cmp"]);
        double result = c->cmp(symtmp,tmp);
        delete tmp;
        delete symtmp;

        //cout << result << endl;
        return result;
}