#include "emfft.h"
#include "cmp.h"
#include "aligner.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 double | refalifnfast (const gsl_vector *v, void *params) |
static Transform | refalin3d_perturb (const Transform *const t, const float &delta, const float &arc, const float &phi, const float &x, const float &y, const float &z) |
static double | refalifn3d (const gsl_vector *v, void *params) |
void | EMAN::dump_aligners () |
map< string, vector< string > > | EMAN::dump_aligners_list () |
|
Definition at line 49 of file aligner.cpp. |
|
Definition at line 950 of file aligner.cpp. References EMAN::Cmp::cmp(), EMAN::EMData::process(), EMAN::Transform::set_mirror(), EMAN::Transform::set_trans(), t, x, and y. Referenced by EMAN::RefineAligner::align(). 00951 { 00952 Dict *dict = (Dict *) params; 00953 00954 double x = gsl_vector_get(v, 0); 00955 double y = gsl_vector_get(v, 1); 00956 double a = gsl_vector_get(v, 2); 00957 00958 EMData *this_img = (*dict)["this"]; 00959 EMData *with = (*dict)["with"]; 00960 bool mirror = (*dict)["mirror"]; 00961 00962 // float mean = (float)this_img->get_attr("mean"); 00963 // if ( Util::goodf(&mean) ) { 00964 // //cout << "tmps mean is nan even before rotation" << endl; 00965 // } 00966 00967 Transform t(Dict("type","2d","alpha",static_cast<float>(a))); 00968 // Transform3D t3d(Transform3D::EMAN, (float)a, 0.0f, 0.0f); 00969 // t3d.set_posttrans( (float) x, (float) y); 00970 // tmp->rotate_translate(t3d); 00971 t.set_trans((float)x,(float)y); 00972 t.set_mirror(mirror); 00973 EMData *tmp = this_img->process("xform",Dict("transform",&t)); 00974 00975 Cmp* c = (Cmp*) ((void*)(*dict)["cmp"]); 00976 double result = c->cmp(tmp,with); 00977 00978 // DELETE AT SOME STAGE, USEFUL FOR PRERELEASE STUFF 00979 // float test_result = (float)result; 00980 // if ( Util::goodf(&test_result) ) { 00981 // cout << "result " << result << " " << x << " " << y << " " << a << endl; 00982 // cout << (float)this_img->get_attr("mean") << " " << (float)tmp->get_attr("mean") << " " << (float)with->get_attr("mean") << endl; 00983 // tmp->write_image("tmp.hdf"); 00984 // with->write_image("with.hdf"); 00985 // this_img->write_image("this_img.hdf"); 00986 // EMData* t = this_img->copy(); 00987 // cout << (float)t->get_attr("mean") << endl; 00988 // t->rotate_translate( t3d ); 00989 // cout << (float)t->get_attr("mean") << endl; 00990 // cout << "exit" << endl; 00992 // cout << (float)t->get_attr("mean") << endl; 00993 // cout << "now exit" << endl; 00994 // delete t; 00995 // } 00996 00997 00998 if ( tmp != 0 ) delete tmp; 00999 01000 return result; 01001 }
|
|
Definition at line 1207 of file aligner.cpp. References EMAN::Cmp::cmp(), phi, EMAN::EMData::process(), refalin3d_perturb(), t, x, and y. Referenced by EMAN::Refine3DAligner::align(). 01208 { 01209 Dict *dict = (Dict *) params; 01210 double x = gsl_vector_get(v, 0); 01211 double y = gsl_vector_get(v, 1); 01212 double z = gsl_vector_get(v, 2); 01213 double arc = gsl_vector_get(v, 3); 01214 double delta = gsl_vector_get(v, 4); 01215 double phi = gsl_vector_get(v, 5); 01216 EMData *this_img = (*dict)["this"]; 01217 EMData *with = (*dict)["with"]; 01218 // bool mirror = (*dict)["mirror"]; 01219 01220 Transform* t = (*dict)["transform"]; 01221 01222 Transform soln = refalin3d_perturb(t,(float)delta,(float)arc,(float)phi,(float)x,(float)y,(float)z); 01223 01224 EMData *tmp = this_img->process("xform",Dict("transform",&soln)); 01225 Cmp* c = (Cmp*) ((void*)(*dict)["cmp"]); 01226 double result = c->cmp(tmp,with); 01227 if ( tmp != 0 ) delete tmp; 01228 delete t; t = 0; 01229 // cout << result << " " << az << " " << alt << " " << phi << " " << x << " " << y << " " << z << endl; 01230 return result; 01231 }
|
|
Definition at line 1003 of file aligner.cpp. References EMAN::EMData::dot_rotate_translate(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), x, and y. Referenced by EMAN::RefineAligner::align(). 01004 { 01005 Dict *dict = (Dict *) params; 01006 EMData *this_img = (*dict)["this"]; 01007 EMData *img_to = (*dict)["with"]; 01008 bool mirror = (*dict)["mirror"]; 01009 01010 double x = gsl_vector_get(v, 0); 01011 double y = gsl_vector_get(v, 1); 01012 double a = gsl_vector_get(v, 2); 01013 01014 double r = this_img->dot_rotate_translate(img_to, (float)x, (float)y, (float)a, mirror); 01015 int nsec = this_img->get_xsize() * this_img->get_ysize(); 01016 double result = 1.0 - r / nsec; 01017 01018 // cout << result << " x " << x << " y " << y << " az " << a << endl; 01019 return result; 01020 }
|
|
Definition at line 1152 of file aligner.cpp. References EMAN::Transform::get_params(), EMAN::Vec3< Type >::normalize(), q, EMAN::Transform::set_trans(), and EMAN::Transform::transpose(). Referenced by EMAN::Refine3DAligner::align(), and refalifn3d(). 01153 { 01154 Dict orig_params = t->get_params("eman"); 01155 float orig_phi = orig_params["phi"]; 01156 float orig_x = orig_params["tx"]; 01157 float orig_y = orig_params["ty"]; 01158 float orig_z = orig_params["tz"]; 01159 orig_params["phi"] = 0; 01160 orig_params["tx"] = 0; 01161 orig_params["ty"] = 0; 01162 orig_params["tz"] = 0; 01163 Transform t_no_phi(orig_params); 01164 01165 Vec3f zz(0,0,1); 01166 01167 Vec3f vv = t_no_phi.transpose()*zz; 01168 Vec3f normal = Vec3f(-vv[2],0,-vv[0]); 01169 01170 normal.normalize(); 01171 01172 Dict d; 01173 d["type"] = "spin"; 01174 d["Omega"] = arc; 01175 d["n1"] = vv[0]; 01176 d["n2"] = vv[1]; 01177 d["n3"] = vv[2]; 01178 01179 Transform q(d); 01180 01181 Vec3f rot_vec = q*normal; 01182 rot_vec.normalize(); 01183 01184 Dict e; 01185 e["type"] = "spin"; 01186 e["Omega"] = delta; 01187 e["n1"] = rot_vec[0]; 01188 e["n2"] = rot_vec[1]; 01189 e["n3"] = rot_vec[2]; 01190 01191 Transform perturb(e); 01192 01193 Dict g; 01194 g["type"] = "eman"; 01195 g["alt"] = 0; 01196 g["az"] = 0; 01197 g["phi"] = 0*phi+orig_phi; 01198 01199 Transform phi_rot(g); 01200 Transform soln = t_no_phi*perturb*phi_rot; 01201 soln.set_trans(x+orig_x,y+orig_y,z+orig_z); 01202 01203 Dict params = soln.get_params("eman"); 01204 return soln; 01205 }
|