#include <random>
#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"
Go to the source code of this file.
|
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) |
|
◆ EMAN2_ALIGNER_DEBUG
#define EMAN2_ALIGNER_DEBUG 0 |
◆ refalidf()
static void refalidf |
( |
const gsl_vector * |
v, |
|
|
void * |
params, |
|
|
gsl_vector * |
df |
|
) |
| |
|
static |
Definition at line 1720 of file aligner.cpp.
1724 static double lstep[4] = { 0.05, 0.05, 0.1, 0.01 };
1726 gsl_vector *vc = gsl_vector_alloc(v->size);
1727 gsl_vector_memcpy(vc,v);
1730 for (
unsigned int i=0; i<v->size; i++) {
1732 double vp = gsl_vector_get(vc,i);
1734 gsl_vector_set(vc,i,
vp+lstep[i]);
1737 gsl_vector_set(vc,i,
vp);
1739 gsl_vector_set(df,i,(f2-f)/lstep[i]);
1742 gsl_vector_free(vc);
static double refalifn(const gsl_vector *v, void *params)
References refalifn().
◆ refalifdf()
static void refalifdf |
( |
const gsl_vector * |
v, |
|
|
void * |
params, |
|
|
double * |
f, |
|
|
gsl_vector * |
df |
|
) |
| |
|
static |
Definition at line 1746 of file aligner.cpp.
1750 static double lstep[4] = { 0.05, 0.05, 0.1, 0.01 };
1752 gsl_vector *vc = gsl_vector_alloc(v->size);
1753 gsl_vector_memcpy(vc,v);
1756 for (
unsigned int i=0; i<v->size; i++) {
1758 double vp = gsl_vector_get(vc,i);
1760 gsl_vector_set(vc,i,
vp+lstep[i]);
1763 gsl_vector_set(vc,i,
vp);
1765 gsl_vector_set(df,i,(f2-*f)/lstep[i]);
1768 gsl_vector_free(vc);
References refalifn().
◆ refalifn()
static double refalifn |
( |
const gsl_vector * |
v, |
|
|
void * |
params |
|
) |
| |
|
static |
Definition at line 1688 of file aligner.cpp.
1692 double x = gsl_vector_get(v, 0);
1693 double y = gsl_vector_get(v, 1);
1694 double a = gsl_vector_get(v, 2);
1696 EMData *this_img = (*dict)[
"this"];
1697 EMData *with = (*dict)[
"with"];
1698 bool mirror = (*dict)[
"mirror"];
1700 Transform t(
Dict(
"type",
"2d",
"alpha",
static_cast<float>(a)));
1701 t.set_trans((
float)
x,(
float)
y);
1702 t.set_mirror(mirror);
1704 float sca=(float)gsl_vector_get(v, 3);
1705 if (sca<.7 || sca>1.3)
return 1.0e20;
1706 t.set_scale((
float)gsl_vector_get(v, 3));
1708 EMData *tmp = this_img->process(
"xform",
Dict(
"transform",&t));
1709 if (dict->
has_key(
"mask")) tmp->mult(*(
EMData *)((*dict)[
"mask"]));
1712 Cmp* c = (
Cmp*) ((
void*)(*dict)[
"cmp"]);
1713 double result = c->
cmp(tmp,with);
1715 if (tmp != 0)
delete tmp;
Cmp class defines image comparison method.
virtual float cmp(EMData *image, EMData *with) const =0
To compare 'image' with another image passed in through its parameters.
Dict is a dictionary to store <string, EMObject> pair.
bool has_key(const string &key) const
Ask the Dictionary if it as a particular key.
EMData stores an image's data and defines core image processing routines.
References EMAN::Cmp::cmp(), EMAN::Dict::has_key(), EMAN::Transform::set_mirror(), EMAN::Transform::set_scale(), EMAN::Transform::set_trans(), x, and y.
Referenced by refalidf(), and refalifdf().
◆ refalifn3dquat()
static double refalifn3dquat |
( |
const gsl_vector * |
v, |
|
|
void * |
params |
|
) |
| |
|
static |
Definition at line 2104 of file aligner.cpp.
2108 double n0 = gsl_vector_get(v, 0);
2109 double n1 = gsl_vector_get(v, 1);
2110 double n2 = gsl_vector_get(v, 2);
2111 double x = gsl_vector_get(v, 3);
2112 double y = gsl_vector_get(v, 4);
2113 double z = gsl_vector_get(v, 5);
2115 EMData *this_img = (*dict)[
"this"];
2116 EMData *with = (*dict)[
"with"];
2119 float spincoeff = (*dict)[
"spincoeff"];
2123 EMData *tmp = this_img->process(
"xform",
Dict(
"transform",&soln));
2124 Cmp* c = (
Cmp*) ((
void*)(*dict)[
"cmp"]);
2125 double result = c->
cmp(tmp,with);
2126 if ( tmp != 0 )
delete tmp;
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)
References EMAN::Cmp::cmp(), refalin3d_perturbquat(), x, and y.
◆ refalifnfast()
static double refalifnfast |
( |
const gsl_vector * |
v, |
|
|
void * |
params |
|
) |
| |
|
static |
Definition at line 1772 of file aligner.cpp.
1775 EMData *this_img = (*dict)[
"this"];
1776 EMData *img_to = (*dict)[
"with"];
1777 bool mirror = (*dict)[
"mirror"];
1779 double x = gsl_vector_get(v, 0);
1780 double y = gsl_vector_get(v, 1);
1781 double a = gsl_vector_get(v, 2);
1784 int nsec = this_img->get_xsize() * this_img->get_ysize();
1785 double result = 1.0 - r / nsec;
double dot_rotate_translate(EMData *with, float dx, float dy, float da, const bool mirror=false)
dot product of 2 images.
References EMAN::EMData::dot_rotate_translate(), x, and y.
◆ refalin3d_perturbquat()
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 |
◆ symquat()
static double symquat |
( |
const gsl_vector * |
v, |
|
|
void * |
params |
|
) |
| |
|
static |
Definition at line 2076 of file aligner.cpp.
2080 double n0 = gsl_vector_get(v, 0);
2081 double n1 = gsl_vector_get(v, 1);
2082 double n2 = gsl_vector_get(v, 2);
2083 double x = gsl_vector_get(v, 3);
2084 double y = gsl_vector_get(v, 4);
2085 double z = gsl_vector_get(v, 5);
2087 EMData* volume = (*dict)[
"volume"];
2088 float spincoeff = (*dict)[
"spincoeff"];
2093 EMData *tmp = volume->process(
"xform",
Dict(
"transform",&soln));
2094 EMData *symtmp = tmp->process(
"xform.applysym",
Dict(
"sym",(*dict)[
"sym"]));
2095 Cmp* c = (
Cmp*) ((
void*)(*dict)[
"cmp"]);
2096 double result = c->
cmp(symtmp,tmp);
References EMAN::Cmp::cmp(), refalin3d_perturbquat(), x, and y.