#include <transform.h>
Public Member Functions | |
Transform () | |
Default constructor Internal matrix is the identity. | |
Transform (const Transform &rhs) | |
Copy constructor. | |
Transform & | operator= (const Transform &that) |
Assignment operator. | |
Transform (const Dict &d) | |
Construction using a dictionary. | |
Transform (const float array[12]) | |
Construction using an array of floats. | |
Transform (const vector< float > array) | |
Construction using a vector of size 12. | |
~Transform () | |
void | set_rotation (const Dict &rotation) |
Set a rotation using a specific Euler type and the dictionary interface Works for all Euler types. | |
void | set_rotation (const Vec3f &v) |
Determine the rotation that would transform a vector pointing in the Z direction so that it points in the direction of the argument vector Automatically normalizes the vector. | |
Dict | get_rotation (const string &euler_type="eman") const |
Get a rotation in any Euler format. | |
Transform | get_rotation_transform () const |
Get the rotation part of the tranformation matrix as a Transform object. | |
Transform | get_hflip_transform () const |
How do I get the transform that will yield the horizontally flipped projection? | |
Transform | get_vflip_transform () const |
How do I get the transform that will yield the vertically flipped projection? | |
void | set_params (const Dict &d) |
Set the parameters of the entire transform. | |
void | set_params_inverse (const Dict &d) |
Set the parameters of the entire transform as though they there in the inverse format. | |
Dict | get_params (const string &euler_type) const |
Get the parameters of the entire transform, using a specific euler convention. | |
Dict | get_params_inverse (const string &euler_type) const |
Get the parameters of the inverse of the transform as though it were in RSMT order not MTSR. | |
void | set_trans (const float &x, const float &y, const float &z=0) |
Set the post translation component. | |
void | set_trans (const Vec3f &v) |
Set the post translation component using a Vec3f. | |
void | set_trans (const Vec2f &v) |
Set the post translation component using a Vec2f. | |
Vec3f | get_trans () const |
Get the post trans as a vec3f. | |
Vec2f | get_trans_2d () const |
Get the degenerant 2D post trans as a vec2f. | |
Vec3f | get_pre_trans () const |
Get the translation vector as though this object was MSRT_ not MTSR, where T_ is what you want Note M means post x mirror, T means translation, S means scale, and R means rotaiton. | |
Vec2f | get_pre_trans_2d () const |
2D version of getting the translation vector as though this object was MSRT_ not MTSR, where T_ is what you want Note M means post x mirror, T means translation, S means scale, and R means rotation | |
template<typename type> | |
void | set_pre_trans (const type &v) |
Set the translational component of the matrix as though it was MSRT_ not MTSR, where T_ is the pre translation. | |
void | set_scale (const float &scale) |
Set the scale. | |
float | get_scale () const |
Get the scale that was applied. | |
bool | get_mirror () const |
Query whether x_mirroring is occuring. | |
void | set_mirror (const bool x_mirror) |
Set whether or not x_mirroring is occuring. | |
void | get_scale_and_mirror (float &scale, bool &x_mirror) const |
Get scale and x_mirror with 1 function call. | |
void | to_identity () |
Force the internal matrix to become the identity. | |
bool | is_identity () const |
Returns whethers or this matrix is the identity. | |
void | orthogonalize () |
Reorthogonalize the rotation part of the matrix in place. | |
float | get_determinant () const |
Get the determinant of the matrix. | |
void | printme () const |
Print the contents of the internal matrix verbatim to standard out. | |
void | set_matrix (const vector< float > &v) |
Set the transformation matrix using a vector. | |
vector< float > | get_matrix () const |
Get the transformation matrix using a vector. | |
void | invert () |
Get the inverse of this transformation matrix. | |
Transform | inverse () const |
Get the inverse of this transformation matrix. | |
void | transpose_inplace () |
Get the transpose of this transformation matrix. | |
Transform | transpose () const |
Get the transpose of this transformation matrix. | |
float | at (int r, int c) const |
Get the value stored in the internal transformation matrix at at coordinate (r,c). | |
void | set (int r, int c, float value) |
Set the value stored in the internal transformation matrix at at coordinate (r,c) to value. | |
float * | operator[] (int i) |
Operator[] convenience so Transform3D[2][2] etc terminology can be used. | |
const float * | operator[] (int i) const |
Operator[] convenience so Transform3D[2][2] etc terminology can be used. | |
Vec2f | transform (const float &x, const float &y) const |
Transform 2D coordinates using the internal transformation matrix. | |
template<typename Type> | |
Vec2f | transform (const Vec2< Type > &v) const |
Transform a 2D vector using the internal transformation matrix. | |
Vec3f | transform (const float &x, const float &y, const float &z) const |
Transform 3D coordinates using the internal transformation matrix. | |
template<typename Type> | |
Vec3f | transform (const Vec3< Type > &v) const |
Transform a 3D vector using the internal transformation matrix. | |
Vec3f | get_matrix3_row (int i) const |
Get a matrix row as a Vec3f required for back compatibility with Tranform3D - see PawelProjector. | |
Transform | get_sym (const string &sym, int n) const |
Apply the symmetry deduced from the function arguments to this Transform and return the result. | |
void | copy_matrix_into_array (float *const) const |
Transform | negate () const |
Negates the Transform - a useful way, for example, for getting an orientation on the opposite side of the sphere. | |
Static Public Member Functions | |
static int | get_nsym (const string &sym) |
get the number of symmetries associated with the given symmetry name | |
static Transform | tet_3_to_2 () |
Get the transform that moves any tetrahedron generated by eman2 so that it matches the 2-2-2 (MRC, FREALIGN) convention. | |
static Transform | icos_5_to_2 () |
Get the transform that moves any icosahedron generated by eman2 so that it matches the 2-2-2 (MRC, FREALIGN) convention. | |
Static Public Attributes | |
static const float | ERR_LIMIT = 0.000001f |
Private Member Functions | |
void | assert_valid_2d () const |
void | init_permissable_keys () |
Called internally to initialize permissable_2d_not_rot, permissable_3d_not_rot, and permissable_rot_keys static members. | |
void | detect_problem_keys (const Dict &d) |
Test to ensure the parametes in the given dictionary are valid Throws if an error is detected Generic - works in every circumstance (set_params, set_rotation, set_params_inv) Uses static members permissable_2d_not_rot, permissable_3d_not_rot, and permissable_rot_keys as basis of decision. | |
Private Attributes | |
float | matrix [3][4] |
Static Private Attributes | |
static vector< string > | permissable_2d_not_rot |
This map is used to validate keys in the argument maps - e.g. if the type is 2d and the angle is not "alpha" then we should throw. | |
static vector< string > | permissable_3d_not_rot |
static map< string, vector< string > > | permissable_rot_keys |
It's designed to store four transformations in a specific order, namely Transform = MTSR Where M is a mirroring operation (about the x-axis) or the identity T is a Translation matrix S is a uniform scaling matrix R is a rotation matrix This means you can call set_scale, set_trans, set_rotation in any order but still have the operations arranged internally in the order of MTSR. This is somewhat restrictive, for example in the context of how OpenGL handles transformations, but in practice is nicely suited to the situations that arise in EMAN2 - namely, alignment and projection orientation characterization.
Note that you can fool the Transform object into storing any matrix by using the constructors that take array arguments. This can useful, for example, for shearing your image.
See http://blake.bcm.tmc.edu/emanwiki/Eman2TransformInPython for using it from Python and detailed discussion See test_transform.py for examples of the way it is unit tested See http://blake.bcm.tmc.edu/emanwiki/EMAN2/Tutorials/RotateTranslate for examples showing how to transform EMDatas with it.
Definition at line 84 of file transform.h.
Transform::Transform | ( | ) |
Default constructor Internal matrix is the identity.
Definition at line 98 of file transform.cpp.
References to_identity().
Referenced by set_params_inverse().
00099 { 00100 to_identity(); 00101 }
Transform::Transform | ( | const Transform & | rhs | ) |
Copy constructor.
rhs | the object to be copied |
Definition at line 103 of file transform.cpp.
Transform::Transform | ( | const Dict & | d | ) |
Construction using a dictionary.
d | the dictionary containing the parameters |
Definition at line 117 of file transform.cpp.
References set_params(), and to_identity().
00117 { 00118 to_identity(); 00119 set_params(d); 00120 }
Transform::Transform | ( | const float | array[12] | ) |
Construction using an array of floats.
array | the array of values that will become the internal matrix. row order (3 rows of 4) |
Definition at line 123 of file transform.cpp.
References matrix.
00123 { 00124 memcpy(matrix,array,12*sizeof(float)); 00125 }
Transform::Transform | ( | const vector< float > | array | ) |
Construction using a vector of size 12.
array | the array of values that will become the internal matrix. row order (3 rows of 4) |
Definition at line 127 of file transform.cpp.
References set_matrix().
00128 { 00129 set_matrix(array); 00130 }
EMAN::Transform::~Transform | ( | ) | [inline] |
void Transform::assert_valid_2d | ( | ) | const [private] |
Definition at line 1162 of file transform.cpp.
References ERR_LIMIT, matrix, and UnexpectedBehaviorException.
Referenced by get_rotation(), and set_rotation().
01162 { 01163 int rotation_error = 0; 01164 int translation_error = 0; 01165 if (fabs(matrix[2][0]) > ERR_LIMIT) rotation_error++; 01166 if (fabs(matrix[2][1]) > ERR_LIMIT) rotation_error++; 01167 if (fabs(matrix[2][3]) > ERR_LIMIT) translation_error++; 01168 if (fabs(matrix[0][2]) > ERR_LIMIT) rotation_error++; 01169 if (fabs(matrix[1][2]) > ERR_LIMIT) rotation_error++; 01170 if ( translation_error && rotation_error ) { 01171 throw UnexpectedBehaviorException("Error, the internal matrix contains 3D rotations and 3D translations. This object can not be considered 2D"); 01172 } else if ( translation_error ) { 01173 throw UnexpectedBehaviorException("Error, the internal matrix contains a non zero z component for a 3D translation. This object can not be considered 2D"); 01174 } 01175 else if ( rotation_error ) { 01176 throw UnexpectedBehaviorException("Error, the internal matrix contains 3D rotations and this object can not be considered 2D"); 01177 } 01178 01179 }
float EMAN::Transform::at | ( | int | r, | |
int | c | |||
) | const [inline] |
Get the value stored in the internal transformation matrix at at coordinate (r,c).
Definition at line 325 of file transform.h.
References matrix.
Referenced by EMAN::PawelProjector::backproject3d().
00325 { return matrix[r][c]; }
void Transform::copy_matrix_into_array | ( | float * | const | ) | const |
Definition at line 143 of file transform.cpp.
References matrix.
Referenced by EMAN::TransformProcessor::process(), and EMAN::TransformProcessor::process_inplace().
00143 { 00144 00145 int idx = 0; 00146 for(int i=0; i<3; ++i) { 00147 for(int j=0; j<4; ++j) { 00148 array[idx] = matrix[i][j]; 00149 idx ++; 00150 } 00151 } 00152 }
void Transform::detect_problem_keys | ( | const Dict & | d | ) | [private] |
Test to ensure the parametes in the given dictionary are valid Throws if an error is detected Generic - works in every circumstance (set_params, set_rotation, set_params_inv) Uses static members permissable_2d_not_rot, permissable_3d_not_rot, and permissable_rot_keys as basis of decision.
d | the dictionary that was the function argument of the set_params, set_rotation or the set_params_inv function |
InvalidParameterException | if the dictionary is invalid in anyway |
Definition at line 313 of file transform.cpp.
References EMAN::Dict::begin(), copy(), EMAN::Dict::end(), EMAN::Dict::has_key_ci(), init_permissable_keys(), InvalidParameterException, permissable_2d_not_rot, permissable_3d_not_rot, permissable_rot_keys, and EMAN::Util::str_to_lower().
Referenced by set_params(), set_params_inverse(), and set_rotation().
00313 { 00314 if (permissable_rot_keys.size() == 0 ) { 00315 init_permissable_keys(); 00316 } 00317 00318 vector<string> verification; 00319 vector<string> problem_keys; 00320 bool is_2d = false; 00321 if (d.has_key_ci("type") ) { 00322 string type = Util::str_to_lower((string)d["type"]); 00323 bool problem = false; 00324 if (permissable_rot_keys.find(type) == permissable_rot_keys.end() ) { 00325 problem_keys.push_back(type); 00326 problem = true; 00327 } 00328 if ( !problem ) { 00329 vector<string> perm = permissable_rot_keys[type]; 00330 std::copy(perm.begin(),perm.end(),back_inserter(verification)); 00331 00332 if ( type == "2d" ) { 00333 is_2d = true; 00334 std::copy(permissable_2d_not_rot.begin(),permissable_2d_not_rot.end(),back_inserter(verification)); 00335 } 00336 } 00337 } 00338 if ( !is_2d ) { 00339 std::copy(permissable_3d_not_rot.begin(),permissable_3d_not_rot.end(),back_inserter(verification)); 00340 } 00341 00342 for (Dict::const_iterator it = d.begin(); it != d.end(); ++it) { 00343 if ( std::find(verification.begin(),verification.end(), it->first) == verification.end() ) { 00344 problem_keys.push_back(it->first); 00345 } 00346 } 00347 00348 if (problem_keys.size() != 0 ) { 00349 string error; 00350 if (problem_keys.size() == 1) { 00351 error = "Transform Error: The \"" +problem_keys[0]+ "\" key is unsupported"; 00352 } else { 00353 error = "Transform Error: The "; 00354 for(vector<string>::const_iterator cit = problem_keys.begin(); cit != problem_keys.end(); ++cit ) { 00355 if ( cit != problem_keys.begin() ) { 00356 if (cit == (problem_keys.end() -1) ) error += " and "; 00357 else error += ", "; 00358 } 00359 error += "\""; 00360 error += *cit; 00361 error += "\""; 00362 } 00363 error += " keys are unsupported"; 00364 } 00365 throw InvalidParameterException(error); 00366 } 00367 }
float Transform::get_determinant | ( | ) | const |
Get the determinant of the matrix.
Definition at line 1078 of file transform.cpp.
References EMAN::Util::apply_precision(), ERR_LIMIT, and matrix.
Referenced by get_mirror(), get_scale(), and get_scale_and_mirror().
01079 { 01080 float det = matrix[0][0]*(matrix[1][1]*matrix[2][2]-matrix[2][1]*matrix[1][2]); 01081 det -= matrix[0][1]*(matrix[1][0]*matrix[2][2]-matrix[2][0]*matrix[1][2]); 01082 det += matrix[0][2]*(matrix[1][0]*matrix[2][1]-matrix[2][0]*matrix[1][1]); 01083 01084 Util::apply_precision(det,ERR_LIMIT); 01085 01086 return det; 01087 }
Transform Transform::get_hflip_transform | ( | ) | const |
How do I get the transform that will yield the horizontally flipped projection?
Definition at line 671 of file transform.cpp.
References get_rotation(), and get_trans().
00671 { 00672 00673 Dict rot = get_rotation("eman"); 00674 rot["alt"] = 180.0f + static_cast<float>(rot["alt"]); 00675 rot["phi"] = 180.0f - static_cast<float>(rot["phi"]); 00676 00677 00678 00679 Transform ret(*this); // Is the identity 00680 ret.set_rotation(rot); 00681 00682 Vec3f trans = get_trans(); 00683 trans[0] = -trans[0]; 00684 ret.set_trans(trans); 00685 00686 // ret.set_mirror(self.get_mirror()); 00687 00688 return ret; 00689 }
vector< float > Transform::get_matrix | ( | ) | const |
Get the transformation matrix using a vector.
Definition at line 154 of file transform.cpp.
References matrix.
00155 { 00156 vector<float> ret(12); 00157 for(int i=0; i<3; ++i) { 00158 for(int j=0; j<4; ++j) { 00159 ret[i*4+j] = matrix[i][j]; 00160 } 00161 } 00162 return ret; 00163 }
Vec3f EMAN::Transform::get_matrix3_row | ( | int | i | ) | const [inline] |
Get a matrix row as a Vec3f required for back compatibility with Tranform3D - see PawelProjector.
i | the row number (starting at 0) |
Definition at line 394 of file transform.h.
References matrix.
Referenced by EMAN::PawelProjector::project3d().
bool Transform::get_mirror | ( | ) | const |
Query whether x_mirroring is occuring.
Definition at line 1049 of file transform.cpp.
References get_determinant().
Referenced by get_params(), get_params_inverse(), get_trans(), get_trans_2d(), EMAN::BackProjectionReconstructor::preprocess_slice(), EMAN::GaussFFTProjector::project3d(), set_mirror(), and set_trans().
01049 { 01050 float determinant = get_determinant(); 01051 01052 bool x_mirror = false; 01053 if ( determinant < 0 ) x_mirror = true; 01054 01055 return x_mirror; 01056 01057 }
int Transform::get_nsym | ( | const string & | sym | ) | [static] |
get the number of symmetries associated with the given symmetry name
Definition at line 1192 of file transform.cpp.
References EMAN::Symmetry3D::get_nsym().
Referenced by EMAN::nnSSNR_ctfReconstructor::setup(), EMAN::nn4_ctfReconstructor::setup(), EMAN::nnSSNR_Reconstructor::setup(), EMAN::nn4Reconstructor::setup(), and EMAN::EMData::symvol().
01193 { 01194 Symmetry3D* sym = Factory<Symmetry3D>::get(sym_name); 01195 int nsym = sym->get_nsym(); 01196 delete sym; 01197 return nsym; 01198 }
Dict Transform::get_params | ( | const string & | euler_type | ) | const |
Get the parameters of the entire transform, using a specific euler convention.
euler_type | the euler type of the retrieved rotation |
Definition at line 423 of file transform.cpp.
References get_mirror(), get_rotation(), get_scale(), get_trans(), EMAN::Util::str_to_lower(), and v.
Referenced by refalin3d_perturb().
00423 { 00424 Dict params = get_rotation(euler_type); 00425 00426 Vec3f v = get_trans(); 00427 params["tx"] = v[0]; params["ty"] = v[1]; 00428 00429 string type = Util::str_to_lower(euler_type); 00430 if ( type != "2d") params["tz"] = v[2]; 00431 00432 float scale = get_scale(); 00433 params["scale"] = scale; 00434 00435 bool mirror = get_mirror(); 00436 params["mirror"] = mirror; 00437 00438 return params; 00439 }
Dict Transform::get_params_inverse | ( | const string & | euler_type | ) | const |
Get the parameters of the inverse of the transform as though it were in RSMT order not MTSR.
euler_type | the euler type of the retrieved rotation |
Definition at line 443 of file transform.cpp.
References get_mirror(), get_pre_trans(), get_rotation(), get_scale(), inverse(), EMAN::Util::str_to_lower(), and v.
00443 { 00444 Transform inv(inverse()); 00445 00446 Dict params = inv.get_rotation(euler_type); 00447 Vec3f v = inv.get_pre_trans(); 00448 params["tx"] = v[0]; params["ty"] = v[1]; 00449 00450 string type = Util::str_to_lower(euler_type); 00451 if ( type != "2d") params["tz"] = v[2]; 00452 00453 float scale = inv.get_scale(); 00454 params["scale"] = scale; 00455 00456 bool mirror = inv.get_mirror(); 00457 params["mirror"] = mirror; 00458 00459 return params; 00460 }
Vec3f Transform::get_pre_trans | ( | ) | const |
Get the translation vector as though this object was MSRT_ not MTSR, where T_ is what you want Note M means post x mirror, T means translation, S means scale, and R means rotaiton.
Definition at line 901 of file transform.cpp.
References get_trans(), invert(), and set_trans().
Referenced by get_params_inverse().
00902 { 00903 Transform T(*this); 00904 T.set_trans(0,0,0); 00905 T.invert(); 00906 00907 Transform soln = T*(*this); 00908 // soln.printme(); 00909 return soln.get_trans(); 00910 }
Vec2f Transform::get_pre_trans_2d | ( | ) | const |
2D version of getting the translation vector as though this object was MSRT_ not MTSR, where T_ is what you want Note M means post x mirror, T means translation, S means scale, and R means rotation
Definition at line 912 of file transform.cpp.
References get_trans_2d(), invert(), and set_trans().
00913 { 00914 Transform T(*this); 00915 T.set_trans(0,0,0); 00916 T.invert(); 00917 00918 Transform soln = T*(*this); 00919 // soln.printme(); 00920 return soln.get_trans_2d(); 00921 }
Dict Transform::get_rotation | ( | const string & | euler_type = "eman" |
) | const |
Get a rotation in any Euler format.
euler_type | the requested Euler type |
Definition at line 707 of file transform.cpp.
References assert_valid_2d(), ERR_LIMIT, get_scale_and_mirror(), InvalidStringException, matrix, max, phi, EMAN::EMConsts::rad2deg, sqrt(), EMAN::Util::str_to_lower(), and UnexpectedBehaviorException.
Referenced by EMAN::RotateTranslateAligner::align(), EMAN::RotationalAligner::align(), EMAN::PawelProjector::backproject3d(), EMAN::ChaoProjector::backproject3d(), get_hflip_transform(), get_params(), get_params_inverse(), get_vflip_transform(), EMAN::FourierReconstructorSimple2D::insert_slice(), EMAN::ChaoProjector::project3d(), EMAN::FourierGriddingProjector::project3d(), and set_pre_trans().
00708 { 00709 Dict result; 00710 00711 float max = 1 - ERR_LIMIT; 00712 float scale; 00713 bool x_mirror; 00714 get_scale_and_mirror(scale,x_mirror); 00715 if (scale == 0) throw UnexpectedBehaviorException("The determinant of the Transform is 0. This is unexpected."); 00716 float cosalt=matrix[2][2]/scale; 00717 float x_mirror_scale = (x_mirror ? -1.0f : 1.0f); 00718 float inv_scale = 1.0f/scale; 00719 00720 float az=0; 00721 float alt = 0; 00722 float phi=0; 00723 float phiS = 0; // like az (but in SPIDER ZXZ) 00724 float psiS =0; // like phi (but in SPIDER ZYZ) 00725 00726 // get alt, az, phi in EMAN convention 00727 if (cosalt > max) { // that is, alt close to 0 00728 alt = 0; 00729 az=0; 00730 phi = (float)EMConsts::rad2deg*(float)atan2(x_mirror_scale*matrix[0][1], x_mirror_scale*matrix[0][0]); 00731 } 00732 else if (cosalt < -max) { // alt close to pi 00733 alt = 180; 00734 az=0; 00735 phi=360.0f-(float)EMConsts::rad2deg*(float)atan2(x_mirror_scale*matrix[0][1], x_mirror_scale*matrix[0][0]); 00736 } 00737 else { 00738 alt = (float)EMConsts::rad2deg*(float) acos(cosalt); 00739 az = 360.0f+(float)EMConsts::rad2deg*(float)atan2(matrix[2][0], -matrix[2][1]); 00740 phi = 360.0f+(float)EMConsts::rad2deg*(float)atan2(x_mirror_scale*matrix[0][2], matrix[1][2]); 00741 } 00742 // az=fmod(az+180.0f,360.0f)-180.0f; 00743 // phi=fmod(phi+180.0f,360.0f)-180.0f; 00744 00745 if (phi > 0) phi = fmod(phi,360.f); 00746 else phi = 360.f - fmod(fabs(phi),360.f); 00747 if (phi == 360.f) phi = 0.f; 00748 00749 if (az > 0 ) az = fmod(az,360.f); 00750 else az = 360.f - fmod(fabs(az),360.f); 00751 if (az == 360.f) az = 0.f; 00752 00753 // get phiS, psiS ; SPIDER 00754 if (fabs(cosalt) > max) { // that is, alt close to 0 00755 phiS=0; 00756 psiS = phi; 00757 } 00758 else { 00759 phiS = az - 90.0f; 00760 psiS = phi + 90.0f; 00761 } 00762 phiS = fmod((phiS + 360.0f ), 360.0f) ; 00763 psiS = fmod((psiS + 360.0f ), 360.0f) ; 00764 00765 // do some quaternionic stuff here 00766 00767 float nphi = (az-phi)/2.0f; 00768 // The next is also e0 00769 float cosOover2 = (cos((az+phi)*M_PI/360) * cos(alt*M_PI/360)) ; 00770 float sinOover2 = sqrt(1 -cosOover2*cosOover2); 00771 float cosnTheta = sin((az+phi)*M_PI/360) * cos(alt*M_PI/360) / sqrt(1-cosOover2*cosOover2) ; 00772 float sinnTheta = sqrt(1-cosnTheta*cosnTheta); 00773 float n1 = sinnTheta*cos(nphi*M_PI/180); 00774 float n2 = sinnTheta*sin(nphi*M_PI/180); 00775 float n3 = cosnTheta; 00776 float xtilt = 0; 00777 float ytilt = 0; 00778 float ztilt = 0; 00779 00780 00781 if (cosOover2<0) { 00782 cosOover2*=-1; n1 *=-1; n2*=-1; n3*=-1; 00783 } 00784 00785 string type = Util::str_to_lower(euler_type); 00786 00787 result["type"] = type; 00788 if (type == "2d") { 00789 assert_valid_2d(); 00790 result["alpha"] = phi; 00791 } else if (type == "eman") { 00792 // assert_consistent_type(THREED); 00793 result["az"] = az; 00794 result["alt"] = alt; 00795 result["phi"] = phi; 00796 } else if (type == "imagic") { 00797 // assert_consistent_type(THREED); 00798 result["alpha"] = az; 00799 result["beta"] = alt; 00800 result["gamma"] = phi; 00801 } else if (type == "spider") { 00802 // assert_consistent_type(THREED); 00803 result["phi"] = phiS; // The first Euler like az 00804 result["theta"] = alt; 00805 result["psi"] = psiS; 00806 } else if (type == "mrc") { 00807 // assert_consistent_type(THREED); 00808 result["phi"] = phiS; 00809 result["theta"] = alt; 00810 result["omega"] = psiS; 00811 } else if (type == "xyz") { 00812 // assert_consistent_type(THREED); 00813 xtilt = atan2(-sin((M_PI/180.0f)*phiS)*sin((M_PI/180.0f)*alt),cos((M_PI/180.0f)*alt)); 00814 ytilt = asin( cos((M_PI/180.0f)*phiS)*sin((M_PI/180.0f)*alt)); 00815 ztilt = psiS*M_PI/180.0f - atan2(sin(xtilt), cos(xtilt) *sin(ytilt)); 00816 00817 xtilt=fmod(xtilt*180/M_PI+540.0f,360.0f) -180.0f; 00818 ztilt=fmod(ztilt*180/M_PI+540.0f,360.0f) -180.0f; 00819 00820 result["xtilt"] = xtilt; 00821 result["ytilt"] = ytilt*180/M_PI; 00822 result["ztilt"] = ztilt; 00823 } else if (type == "quaternion") { 00824 // assert_consistent_type(THREED); 00825 result["e0"] = cosOover2 ; 00826 result["e1"] = sinOover2 * n1 ; 00827 result["e2"] = sinOover2 * n2; 00828 result["e3"] = sinOover2 * n3; 00829 } else if (type == "spin") { 00830 // assert_consistent_type(THREED); 00831 result["Omega"] =360.0f* acos(cosOover2)/ M_PI ; 00832 result["n1"] = n1; 00833 result["n2"] = n2; 00834 result["n3"] = n3; 00835 } else if (type == "sgirot") { 00836 // assert_consistent_type(THREED); 00837 result["q"] = 360.0f*acos(cosOover2)/M_PI ; 00838 result["n1"] = n1; 00839 result["n2"] = n2; 00840 result["n3"] = n3; 00841 } else if (type == "matrix") { 00842 // assert_consistent_type(THREED); 00843 result["m11"] = x_mirror_scale*matrix[0][0]*inv_scale; 00844 result["m12"] = x_mirror_scale*matrix[0][1]*inv_scale; 00845 result["m13"] = x_mirror_scale*matrix[0][2]*inv_scale; 00846 result["m21"] = matrix[1][0]*inv_scale; 00847 result["m22"] = matrix[1][1]*inv_scale; 00848 result["m23"] = matrix[1][2]*inv_scale; 00849 result["m31"] = matrix[2][0]*inv_scale; 00850 result["m32"] = matrix[2][1]*inv_scale; 00851 result["m33"] = matrix[2][2]*inv_scale; 00852 } else { 00853 throw InvalidStringException(euler_type, "unknown Euler Type"); 00854 } 00855 00856 return result; 00857 }
Transform Transform::get_rotation_transform | ( | ) | const |
Get the rotation part of the tranformation matrix as a Transform object.
Definition at line 627 of file transform.cpp.
References set_mirror(), set_scale(), and set_trans().
Referenced by EMAN::GaussFFTProjector::project3d().
00628 { 00629 Transform ret(*this); 00630 ret.set_scale(1.0); 00631 ret.set_mirror(false); 00632 ret.set_trans(0,0,0); 00633 //ret.orthogonalize(); // ? 00634 return ret; 00635 }
float Transform::get_scale | ( | ) | const |
Get the scale that was applied.
Definition at line 952 of file transform.cpp.
References EMAN::Util::apply_precision(), ERR_LIMIT, and get_determinant().
Referenced by get_params(), get_params_inverse(), EMAN::BackProjectionReconstructor::preprocess_slice(), EMAN::GaussFFTProjector::project3d(), set_pre_trans(), and set_scale().
00952 { 00953 float determinant = get_determinant(); 00954 if (determinant < 0 ) determinant *= -1; 00955 00956 float scale = std::pow(determinant,1.0f/3.0f); 00957 int int_scale = static_cast<int>(scale); 00958 float scale_residual = scale-static_cast<float>(int_scale); 00959 if ( scale_residual < ERR_LIMIT ) { scale = static_cast<float>(int_scale); }; 00960 00961 Util::apply_precision(scale, ERR_LIMIT); 00962 00963 return scale; 00964 }
void Transform::get_scale_and_mirror | ( | float & | scale, | |
bool & | x_mirror | |||
) | const |
Get scale and x_mirror with 1 function call.
More efficient than calling get_scale and get_x_mirror separately
scale | a reference to the value that will be assigned the scale value | |
x_mirror | a reference to the value that will be assigned the x_mirror value |
Definition at line 1059 of file transform.cpp.
References EMAN::Util::apply_precision(), ERR_LIMIT, and get_determinant().
Referenced by get_rotation(), orthogonalize(), and set_rotation().
01059 { 01060 01061 float determinant = get_determinant(); 01062 x_mirror = false; 01063 if ( determinant < 0 ) { 01064 x_mirror = true; 01065 determinant *= -1; 01066 } 01067 if (determinant != 1 ) { 01068 scale = std::pow(determinant,1.0f/3.0f); 01069 int int_scale = static_cast<int>(scale); 01070 float scale_residual = scale-static_cast<float>(int_scale); 01071 if ( scale_residual < ERR_LIMIT ) { scale = static_cast<float>(int_scale); }; 01072 } 01073 else scale = 1; 01074 01075 Util::apply_precision(scale,ERR_LIMIT); 01076 }
Transform Transform::get_sym | ( | const string & | sym, | |
int | n | |||
) | const |
Apply the symmetry deduced from the function arguments to this Transform and return the result.
Definition at line 1183 of file transform.cpp.
References EMAN::Symmetry3D::get_sym().
Referenced by EMAN::EMData::symvol().
01184 { 01185 Symmetry3D* sym = Factory<Symmetry3D>::get(sym_name); 01186 Transform ret; 01187 ret = (*this) * sym->get_sym(n); 01188 delete sym; 01189 return ret; 01190 }
Vec3f Transform::get_trans | ( | ) | const |
Get the post trans as a vec3f.
Definition at line 872 of file transform.cpp.
References EMAN::Util::apply_precision(), ERR_LIMIT, get_mirror(), matrix, and v.
Referenced by get_hflip_transform(), get_params(), get_pre_trans(), get_vflip_transform(), EMAN::GaussFFTProjector::project3d(), EMAN::EMData::rot_scale_trans(), EMAN::EMData::rot_scale_trans_background(), and set_params_inverse().
00873 { 00874 // No type asserted 00875 bool x_mirror = get_mirror(); 00876 Vec3f v; 00877 if (x_mirror) v[0] = -matrix[0][3]; 00878 else v[0] = matrix[0][3]; 00879 v[1] = matrix[1][3]; 00880 v[2] = matrix[2][3]; 00881 00882 Util::apply_precision(v[0],ERR_LIMIT); 00883 Util::apply_precision(v[1],ERR_LIMIT); 00884 Util::apply_precision(v[2],ERR_LIMIT); 00885 00886 return v; 00887 }
Vec2f Transform::get_trans_2d | ( | ) | const |
Get the degenerant 2D post trans as a vec2f.
Definition at line 889 of file transform.cpp.
References get_mirror(), matrix, and v.
Referenced by get_pre_trans_2d(), and EMAN::BackProjectionReconstructor::preprocess_slice().
00890 { 00891 bool x_mirror = get_mirror(); 00892 Vec2f v; 00893 if (x_mirror) v[0] = -matrix[0][3]; 00894 else v[0] = matrix[0][3]; 00895 v[1] = matrix[1][3]; 00896 return v; 00897 }
Transform Transform::get_vflip_transform | ( | ) | const |
How do I get the transform that will yield the vertically flipped projection?
Definition at line 691 of file transform.cpp.
References get_rotation(), and get_trans().
00691 { 00692 00693 Dict rot = get_rotation("eman"); 00694 rot["alt"] = 180.0f + static_cast<float>(rot["alt"]); 00695 rot["phi"] = - static_cast<float>(rot["phi"]); 00696 00697 Transform ret(*this); 00698 ret.set_rotation(rot); 00699 00700 Vec3f trans = get_trans(); 00701 trans[1] = -trans[1]; 00702 ret.set_trans(trans); 00703 00704 return ret; 00705 }
Transform Transform::icos_5_to_2 | ( | ) | [static] |
Get the transform that moves any icosahedron generated by eman2 so that it matches the 2-2-2 (MRC, FREALIGN) convention.
Doctor Phil says: alt = (acos(cos(pi/5)/sqrt(3)/sin(pi/5)) + acos(2*cos(pi/5)/ sqrt(3) ) )*180/pi This is the angle between a 5 and a 3 plus the angle between a 3 and a 2
Definition at line 67 of file transform.cpp.
References t.
00067 { 00068 Transform t; 00069 Dict d; 00070 d["type"] = "eman"; 00071 d["phi"] = 0.0f; 00072 d["az"] = 270.0f; 00073 d["alt"] = 58.282523f; // 5 fold to a 3 fold to a 2 fold 00079 t.set_rotation(d); 00080 return t; 00081 }
void Transform::init_permissable_keys | ( | ) | [private] |
Called internally to initialize permissable_2d_not_rot, permissable_3d_not_rot, and permissable_rot_keys static members.
Definition at line 228 of file transform.cpp.
References permissable_2d_not_rot, permissable_3d_not_rot, and permissable_rot_keys.
Referenced by detect_problem_keys().
00229 { 00230 00231 permissable_2d_not_rot.push_back("tx"); 00232 permissable_2d_not_rot.push_back("ty"); 00233 permissable_2d_not_rot.push_back("scale"); 00234 permissable_2d_not_rot.push_back("mirror"); 00235 permissable_2d_not_rot.push_back("type"); 00236 00237 permissable_3d_not_rot.push_back("tx"); 00238 permissable_3d_not_rot.push_back("ty"); 00239 permissable_3d_not_rot.push_back("tz"); 00240 permissable_3d_not_rot.push_back("scale"); 00241 permissable_3d_not_rot.push_back("mirror"); 00242 permissable_3d_not_rot.push_back("type"); 00243 00244 vector<string> tmp; 00245 tmp.push_back("alpha"); 00246 permissable_rot_keys["2d"] = tmp; 00247 00248 tmp.clear(); 00249 tmp.push_back("alt"); 00250 tmp.push_back("az"); 00251 tmp.push_back("phi"); 00252 permissable_rot_keys["eman"] = tmp; 00253 00254 tmp.clear(); 00255 tmp.push_back("psi"); 00256 tmp.push_back("theta"); 00257 tmp.push_back("phi"); 00258 permissable_rot_keys["spider"] = tmp; 00259 00260 tmp.clear(); 00261 tmp.push_back("alpha"); 00262 tmp.push_back("beta"); 00263 tmp.push_back("gamma"); 00264 permissable_rot_keys["imagic"] = tmp; 00265 00266 tmp.clear(); 00267 tmp.push_back("ztilt"); 00268 tmp.push_back("xtilt"); 00269 tmp.push_back("ytilt"); 00270 permissable_rot_keys["xyz"] = tmp; 00271 00272 tmp.clear(); 00273 tmp.push_back("phi"); 00274 tmp.push_back("theta"); 00275 tmp.push_back("omega"); 00276 permissable_rot_keys["mrc"] = tmp; 00277 00278 tmp.clear(); 00279 tmp.push_back("e0"); 00280 tmp.push_back("e1"); 00281 tmp.push_back("e2"); 00282 tmp.push_back("e3"); 00283 permissable_rot_keys["quaternion"] = tmp; 00284 00285 tmp.clear(); 00286 tmp.push_back("n1"); 00287 tmp.push_back("n2"); 00288 tmp.push_back("n3"); 00289 tmp.push_back("Omega"); 00290 permissable_rot_keys["spin"] = tmp; 00291 00292 tmp.clear(); 00293 tmp.push_back("n1"); 00294 tmp.push_back("n2"); 00295 tmp.push_back("n3"); 00296 tmp.push_back("q"); 00297 permissable_rot_keys["sgirot"] = tmp; 00298 00299 tmp.clear(); 00300 tmp.push_back("m11"); 00301 tmp.push_back("m12"); 00302 tmp.push_back("m13"); 00303 tmp.push_back("m21"); 00304 tmp.push_back("m22"); 00305 tmp.push_back("m23"); 00306 tmp.push_back("m31"); 00307 tmp.push_back("m32"); 00308 tmp.push_back("m33"); 00309 permissable_rot_keys["matrix"] = tmp; 00310 }
Transform Transform::inverse | ( | ) | const |
Get the inverse of this transformation matrix.
Definition at line 1123 of file transform.cpp.
References t.
Referenced by get_params_inverse(), EMAN::EMData::max_3D_pixel_error(), EMAN::StandardProjector::project3d(), EMAN::EMData::rot_scale_trans(), and EMAN::EMData::rot_scale_trans_background().
void Transform::invert | ( | ) |
Get the inverse of this transformation matrix.
Definition at line 1089 of file transform.cpp.
Referenced by EMAN::EMData::extract_plane(), get_pre_trans(), get_pre_trans_2d(), EMAN::Symmetry3D::in_which_asym_unit(), EMAN::BackProjectionReconstructor::insert_slice(), EMAN::GaussFFTProjector::project3d(), EMAN::Symmetry3D::reduce(), set_params_inverse(), and set_pre_trans().
01089 { 01090 01091 float m00 = matrix[0][0]; float m01=matrix[0][1]; float m02=matrix[0][2]; 01092 float m10 = matrix[1][0]; float m11=matrix[1][1]; float m12=matrix[1][2]; 01093 float m20 = matrix[2][0]; float m21=matrix[2][1]; float m22=matrix[2][2]; 01094 float v0 = matrix[0][3]; float v1 =matrix[1][3]; float v2 =matrix[2][3]; 01095 01096 float cof00 = m11*m22-m12*m21; 01097 float cof11 = m22*m00-m20*m02; 01098 float cof22 = m00*m11-m01*m10; 01099 float cof01 = m10*m22-m20*m12; 01100 float cof02 = m10*m21-m20*m11; 01101 float cof12 = m00*m21-m01*m20; 01102 float cof10 = m01*m22-m02*m21; 01103 float cof20 = m01*m12-m02*m11; 01104 float cof21 = m00*m12-m10*m02; 01105 01106 float det = m00* cof00 + m02* cof02 -m01*cof01; 01107 01108 matrix[0][0] = cof00/det; 01109 matrix[0][1] = - cof10/det; 01110 matrix[0][2] = cof20/det; 01111 matrix[1][0] = - cof01/det; 01112 matrix[1][1] = cof11/det; 01113 matrix[1][2] = - cof21/det; 01114 matrix[2][0] = cof02/det; 01115 matrix[2][1] = - cof12/det; 01116 matrix[2][2] = cof22/det; 01117 01118 matrix[0][3] = (- cof00*v0 + cof10*v1 - cof20*v2 )/det; 01119 matrix[1][3] = ( cof01*v0 - cof11*v1 + cof21*v2 )/det; 01120 matrix[2][3] = (- cof02*v0 + cof12*v1 - cof22*v2 )/det; 01121 }
bool Transform::is_identity | ( | ) | const |
Returns whethers or this matrix is the identity.
Definition at line 180 of file transform.cpp.
References EMAN::Util::apply_precision(), ERR_LIMIT, and matrix.
Referenced by EMAN::TestTomoImage::insert_rectangle(), and EMAN::FourierReconstructor::preprocess_slice().
00180 { 00181 for(int i=0; i<3; ++i) { 00182 for(int j=0; j<4; ++j) { 00183 float c = matrix[i][j]; 00184 Util::apply_precision(c,ERR_LIMIT); 00185 if(i==j) { 00186 if (c != 1.0) return false; 00187 } 00188 else { 00189 if (c != 0.0) return false; 00190 } 00191 } 00192 } 00193 return true; 00194 }
Transform Transform::negate | ( | ) | const |
Negates the Transform - a useful way, for example, for getting an orientation on the opposite side of the sphere.
Definition at line 660 of file transform.cpp.
References t.
00661 { 00662 Transform t(*this); 00663 for(unsigned int i = 0; i < 3; ++i) { 00664 for(unsigned int j = 0; j < 4; ++j) { 00665 t.set(i,j,t[i][j]*-1); 00666 } 00667 } 00668 return t; 00669 }
Assignment operator.
that | that which this will become |
Definition at line 108 of file transform.cpp.
References matrix.
00108 { 00109 00110 if (this != &that ) { 00111 memcpy(matrix,that.matrix,12*sizeof(float)); 00112 // transform_type = that.transform_type; 00113 } 00114 return *this; 00115 }
const float* EMAN::Transform::operator[] | ( | int | i | ) | const [inline] |
Operator[] convenience so Transform3D[2][2] etc terminology can be used.
Definition at line 339 of file transform.h.
References matrix.
00339 { return matrix[i]; }
float* EMAN::Transform::operator[] | ( | int | i | ) | [inline] |
Operator[] convenience so Transform3D[2][2] etc terminology can be used.
Definition at line 334 of file transform.h.
References matrix.
00334 { return matrix[i]; }
void Transform::orthogonalize | ( | ) |
Reorthogonalize the rotation part of the matrix in place.
Does this by performing the SVD decomposition of the rotation matrix R such that R = USV^T - since the eigenvalues of a rotation matrix are all 1 we enforce that S should be the identity and produce a corrected matrix R' = UV^T
Definition at line 978 of file transform.cpp.
References get_scale_and_mirror(), matrix, R, UnexpectedBehaviorException, and V.
00979 { 00980 float scale; 00981 bool x_mirror; 00982 get_scale_and_mirror(scale,x_mirror); 00983 if (scale == 0) throw UnexpectedBehaviorException("The determinant of the Transform is 0. This is unexpected."); 00984 double inv_scale = 1.0/static_cast<double>(scale); 00985 double mirror_scale = (x_mirror == true ? -1.0:1.0); 00986 00987 gsl_matrix * R = gsl_matrix_calloc(3,3); 00988 for ( unsigned int i = 0; i < 3; ++i ) 00989 { 00990 for ( unsigned int j = 0; j < 3; ++j ) 00991 { 00992 if (i == 0 && mirror_scale != 1.0 ) { 00993 gsl_matrix_set( R, i, j, static_cast<double>(matrix[i][j])*mirror_scale*inv_scale ); 00994 } 00995 else { 00996 gsl_matrix_set( R, i, j, static_cast<double>(matrix[i][j])*inv_scale ); 00997 } 00998 } 00999 } 01000 01001 gsl_matrix * V = gsl_matrix_calloc(3,3); 01002 gsl_vector * S = gsl_vector_calloc(3); 01003 gsl_vector * work = gsl_vector_calloc(3); 01004 gsl_linalg_SV_decomp (R, V, S, work); // Now R is U of the SVD R = USV^T 01005 01006 gsl_matrix * Soln = gsl_matrix_calloc(3,3); 01007 gsl_blas_dgemm (CblasNoTrans, CblasTrans, 1.0, R, V, 0.0, Soln); 01008 01009 for ( unsigned int i = 0; i < 3; ++i ) 01010 { 01011 for ( unsigned int j = 0; j < 3; ++j ) 01012 { 01013 matrix[i][j] = static_cast<float>( gsl_matrix_get(Soln,i,j) ); 01014 } 01015 } 01016 01017 // Apply scale if it existed previously 01018 if (scale != 1.0f) { 01019 for(int i=0; i<3; ++i) { 01020 for(int j=0; j<3; ++j) { 01021 matrix[i][j] *= scale; 01022 } 01023 } 01024 } 01025 01026 // Apply post x mirroring if it was applied previouslys 01027 if ( x_mirror ) { 01028 for(int j=0; j<3; ++j) { 01029 matrix[0][j] *= -1.0f; 01030 } 01031 } 01032 01033 gsl_matrix_free(V); gsl_matrix_free(R); gsl_matrix_free(Soln); 01034 gsl_vector_free(S); gsl_vector_free(work); 01035 }
void EMAN::Transform::printme | ( | ) | const [inline] |
Print the contents of the internal matrix verbatim to standard out.
Definition at line 286 of file transform.h.
References matrix.
00286 { 00287 cout << matrix[0][0] << " " << matrix[0][1] << " " << matrix[0][2] << " " << matrix[0][3] << endl; 00288 cout << matrix[1][0] << " " << matrix[1][1] << " " << matrix[1][2] << " " << matrix[1][3] << endl; 00289 cout << matrix[2][0] << " " << matrix[2][1] << " " << matrix[2][2] << " " << matrix[2][3] << endl; 00290 cout << 0 << " " << 0 << " " << 0 << " " << 1 << endl; 00291 }
void EMAN::Transform::set | ( | int | r, | |
int | c, | |||
float | value | |||
) | [inline] |
Set the value stored in the internal transformation matrix at at coordinate (r,c) to value.
Definition at line 329 of file transform.h.
References matrix.
00329 { matrix[r][c] = value; }
void Transform::set_matrix | ( | const vector< float > & | v | ) |
Set the transformation matrix using a vector.
Must be of length 12.
v | the transformation matrix stored as a vector - 3 rows of 4. |
Definition at line 132 of file transform.cpp.
References InvalidParameterException, and matrix.
Referenced by EMAN::EMObject::operator Transform *(), and Transform().
00133 { 00134 if (v.size() != 12 ) throw InvalidParameterException("The construction array must be of size 12"); 00135 00136 for(int i=0; i<3; ++i) { 00137 for(int j=0; j<4; ++j) { 00138 matrix[i][j] = v[i*4+j]; 00139 } 00140 } 00141 }
void Transform::set_mirror | ( | const bool | x_mirror | ) |
Set whether or not x_mirroring is occuring.
x_mirror | whether x_mirroring should be applied |
Definition at line 1037 of file transform.cpp.
References get_mirror(), and matrix.
Referenced by EMAN::RefineAligner::align(), EMAN::FourierReconstructor::determine_slice_agreement(), get_rotation_transform(), EMAN::BackProjectionReconstructor::insert_slice(), EMAN::FourierReconstructor::insert_slice(), set_params(), and set_params_inverse().
01037 { 01038 01039 bool old_x_mirror = get_mirror(); 01040 if (old_x_mirror == x_mirror) return; // The user is setting the same value 01041 else { 01042 // Toggle the mirroring operation 01043 for (int j = 0; j < 4; ++j ) { 01044 matrix[0][j] *= -1; 01045 } 01046 } 01047 }
void Transform::set_params | ( | const Dict & | d | ) |
Set the parameters of the entire transform.
keys acted upon are "type" - if this exists then the correct euler angles need to be included - also "tx","ty","tz", "scale", and "mirror"
d | the dictionary containing the parameters |
Definition at line 197 of file transform.cpp.
References EMAN::EMObject::BOOL, detect_problem_keys(), EMAN::Dict::get_ci(), EMAN::Dict::has_key_ci(), EMAN::EMObject::INT, InvalidParameterException, set_mirror(), set_rotation(), set_scale(), set_trans(), and EMAN::EMObject::UNSIGNEDINT.
Referenced by Transform().
00197 { 00198 detect_problem_keys(d); 00199 00200 if (d.has_key_ci("type") ) set_rotation(d); 00201 00202 if (d.has_key_ci("scale")) { 00203 float scale = static_cast<float>(d.get_ci("scale")); 00204 set_scale(scale); 00205 } 00206 00207 float dx=0,dy=0,dz=0; 00208 00209 if (d.has_key_ci("tx")) dx = static_cast<float>(d.get_ci("tx")); 00210 if (d.has_key_ci("ty")) dy = static_cast<float>(d.get_ci("ty")); 00211 if (d.has_key_ci("tz")) dz = static_cast<float>(d.get_ci("tz")); 00212 00213 if ( dx != 0.0 || dy != 0.0 || dz != 0.0 ) { 00214 set_trans(dx,dy,dz); 00215 } 00216 00217 if (d.has_key_ci("mirror")) { 00218 EMObject e = d.get_ci("mirror"); 00219 if ( (e.get_type() != EMObject::BOOL ) && (e.get_type() != EMObject::INT ) && (e.get_type() != EMObject::UNSIGNEDINT ) ) 00220 throw InvalidParameterException("Error, mirror must be a bool or an int"); 00221 00222 bool mirror = static_cast<bool>(e); 00223 set_mirror(mirror); 00224 } 00225 }
void Transform::set_params_inverse | ( | const Dict & | d | ) |
Set the parameters of the entire transform as though they there in the inverse format.
in other words, calling set_params_inverse(get_params_inverse()) should essentially leave the object unchanged.
d | the dictionary containing the inverse parameters |
Definition at line 369 of file transform.cpp.
References EMAN::EMObject::BOOL, detect_problem_keys(), EMAN::Dict::get_ci(), get_trans(), EMAN::Dict::has_key_ci(), EMAN::EMObject::INT, InvalidParameterException, invert(), set_mirror(), set_rotation(), set_scale(), set_trans(), Transform(), and EMAN::EMObject::UNSIGNEDINT.
00369 { 00370 detect_problem_keys(d); 00371 00372 if (d.has_key_ci("type") ) set_rotation(d); 00373 00374 float dx=0,dy=0,dz=0; 00375 if (d.has_key_ci("tx")) dx = static_cast<float>(d.get_ci("tx")); 00376 if (d.has_key_ci("ty")) dy = static_cast<float>(d.get_ci("ty")); 00377 if (d.has_key_ci("tz")) dz = static_cast<float>(d.get_ci("tz")); 00378 00379 if ( (dx != 0.0 || dy != 0.0 || dz != 0.0) && d.has_key_ci("type") ) { 00380 Transform pre_trans; 00381 pre_trans.set_trans(dx,dy,dz); 00382 00383 Transform tmp; 00384 tmp.set_rotation(d); 00385 00386 if (d.has_key_ci("scale")) { 00387 float scale = static_cast<float>(d.get_ci("scale")); 00388 tmp.set_scale(scale); 00389 } 00390 00391 Transform solution_trans = tmp*pre_trans; 00392 00393 if (d.has_key_ci("scale")) { 00394 Transform tmp; 00395 float scale = static_cast<float>(d.get_ci("scale")); 00396 tmp.set_scale(scale); 00397 solution_trans = solution_trans*tmp; 00398 } 00399 00400 tmp = Transform(); 00401 tmp.set_rotation(d); 00402 solution_trans = solution_trans*tmp; 00403 set_trans(solution_trans.get_trans()); 00404 } 00405 00406 if (d.has_key_ci("scale")) { 00407 float scale = static_cast<float>(d.get_ci("scale")); 00408 set_scale(scale); 00409 } 00410 00411 if (d.has_key_ci("mirror")) { 00412 EMObject e = d.get_ci("mirror"); 00413 if ( (e.get_type() != EMObject::BOOL ) && (e.get_type() != EMObject::INT ) && (e.get_type() != EMObject::UNSIGNEDINT ) ) 00414 throw InvalidParameterException("Error, mirror must be a bool or an int"); 00415 00416 bool mirror = static_cast<bool>(e); 00417 set_mirror(mirror); 00418 } 00419 invert(); 00420 }
void EMAN::Transform::set_pre_trans | ( | const type & | v | ) |
Set the translational component of the matrix as though it was MSRT_ not MTSR, where T_ is the pre translation.
Internally the correct form of MTSR is computed.
v | the vector (Vec3f or Vec2f) that is the pre trans |
Definition at line 482 of file transform.h.
References get_rotation(), get_scale(), invert(), set_rotation(), set_scale(), and set_trans().
00482 { 00483 00484 Transform tmp; 00485 Dict rot = get_rotation("eman"); 00486 tmp.set_rotation(rot); 00487 00488 float scale = get_scale(); 00489 if (scale != 1.0 ) tmp.set_scale(scale); 00490 00491 Transform trans; 00492 trans.set_trans(v); 00493 00494 trans = tmp*trans; 00495 00496 Transform tmp2; 00497 tmp2.set_rotation(rot); 00498 tmp2.invert(); // invert 00499 if (scale != 1.0 ) tmp2.set_scale(1.0f/scale); 00500 00501 00502 trans = trans*tmp2; 00503 00504 set_trans(trans.get_trans()); 00505 }
void Transform::set_rotation | ( | const Vec3f & | v | ) |
Determine the rotation that would transform a vector pointing in the Z direction so that it points in the direction of the argument vector Automatically normalizes the vector.
v | the direction you want to solve for |
Definition at line 637 of file transform.cpp.
References EMAN::Vec3< Type >::normalize(), EMAN::EMConsts::rad2deg, set_rotation(), theta, UnexpectedBehaviorException, and v.
00638 { 00639 if ( v[0] == 0 && v[1] == 0 && v[2] == 0 ) 00640 throw UnexpectedBehaviorException("Can't set rotation for the null vector"); 00641 00642 Vec3f v1(v); 00643 v1.normalize(); 00644 00645 float theta = acos(v1[2]); // in radians 00646 float psi = atan2(v1[1],-v1[0]); 00647 // if (v1[1] != 0) psi = asin(v1[1]/sin(theta)); // in radians 00648 00649 Dict d; 00650 d["theta"] = (float)EMConsts::rad2deg*theta; 00651 d["psi"] = (float)EMConsts::rad2deg*psi; 00652 d["phi"] = (float)0.0; 00653 d["type"] = "spider"; 00654 00655 set_rotation(d); 00656 00657 00658 }
void Transform::set_rotation | ( | const Dict & | rotation | ) |
Set a rotation using a specific Euler type and the dictionary interface Works for all Euler types.
rotation | a dictionary containing all key-entry pair required of the associated Euler type |
Definition at line 478 of file transform.cpp.
References assert_valid_2d(), detect_problem_keys(), EMAN::Dict::get_ci(), get_scale_and_mirror(), EMAN::Dict::has_key_ci(), InvalidParameterException, InvalidStringException, matrix, phi, EMAN::Util::str_to_lower(), and UnexpectedBehaviorException.
Referenced by EMAN::RotatePrecenterAligner::align(), EMAN::BackProjectionReconstructor::preprocess_slice(), EMAN::FourierReconstructor::preprocess_slice(), set_params(), set_params_inverse(), set_pre_trans(), and set_rotation().
00479 { 00480 detect_problem_keys(rotation); 00481 string euler_type; 00482 00483 if (!rotation.has_key_ci("type") ){ 00484 throw InvalidParameterException("argument dictionary does not contain the type key"); 00485 } 00486 00487 euler_type = static_cast<string>(rotation.get_ci("type"));// Warning, will throw 00488 00489 00490 float e0=0;float e1=0; float e2=0; float e3=0; 00491 float Omega=0; 00492 float az = 0; 00493 float alt = 0; 00494 float phi = 0; 00495 float cxtilt = 0; 00496 float sxtilt = 0; 00497 float cytilt = 0; 00498 float sytilt = 0; 00499 bool is_quaternion = 0; 00500 bool is_matrix = 0; 00501 00502 bool x_mirror; 00503 float scale; 00504 // Get these before anything changes so we can apply them again after the rotation is set 00505 get_scale_and_mirror(scale,x_mirror); 00506 if (scale == 0) throw UnexpectedBehaviorException("The determinant of the Transform is 0. This is unexpected."); 00507 00508 string type = Util::str_to_lower(euler_type); 00509 if (type == "2d") { 00510 assert_valid_2d(); 00511 az = 0; 00512 alt = 0; 00513 phi = (float)rotation["alpha"] ; 00514 } else if ( type == "eman" ) { 00515 // validate_and_set_type(THREED); 00516 az = (float)rotation["az"] ; 00517 alt = (float)rotation["alt"] ; 00518 phi = (float)rotation["phi"] ; 00519 } else if ( type == "imagic" ) { 00520 // validate_and_set_type(THREED); 00521 az = (float)rotation["alpha"] ; 00522 alt = (float)rotation["beta"] ; 00523 phi = (float)rotation["gamma"] ; 00524 } else if ( type == "spider" ) { 00525 // validate_and_set_type(THREED); 00526 az = (float)rotation["phi"] + 90.0f; 00527 alt = (float)rotation["theta"] ; 00528 phi = (float)rotation["psi"] - 90.0f; 00529 } else if ( type == "xyz" ) { 00530 // validate_and_set_type(THREED); 00531 cxtilt = cos( (M_PI/180.0f)*(float)rotation["xtilt"]); 00532 sxtilt = sin( (M_PI/180.0f)*(float)rotation["xtilt"]); 00533 cytilt = cos( (M_PI/180.0f)*(float)rotation["ytilt"]); 00534 sytilt = sin( (M_PI/180.0f)*(float)rotation["ytilt"]); 00535 az = (180.0f/M_PI)*atan2(-cytilt*sxtilt,sytilt) + 90.0f ; 00536 alt = (180.0f/M_PI)*acos(cytilt*cxtilt) ; 00537 phi = (float)rotation["ztilt"] +(180.0f/M_PI)*atan2(sxtilt,cxtilt*sytilt) - 90.0f ; 00538 } else if ( type == "mrc" ) { 00539 // validate_and_set_type(THREED); 00540 az = (float)rotation["phi"] + 90.0f ; 00541 alt = (float)rotation["theta"] ; 00542 phi = (float)rotation["omega"] - 90.0f ; 00543 } else if ( type == "quaternion" ) { 00544 // validate_and_set_type(THREED); 00545 is_quaternion = 1; 00546 e0 = (float)rotation["e0"]; 00547 e1 = (float)rotation["e1"]; 00548 e2 = (float)rotation["e2"]; 00549 e3 = (float)rotation["e3"]; 00550 } else if ( type == "spin" ) { 00551 // validate_and_set_type(THREED); 00552 is_quaternion = 1; 00553 Omega = (float)rotation["Omega"]; 00554 e0 = cos(Omega*M_PI/360.0f); 00555 e1 = sin(Omega*M_PI/360.0f)* (float)rotation["n1"]; 00556 e2 = sin(Omega*M_PI/360.0f)* (float)rotation["n2"]; 00557 e3 = sin(Omega*M_PI/360.0f)* (float)rotation["n3"]; 00558 } else if ( type == "sgirot" ) { 00559 // validate_and_set_type(THREED); 00560 is_quaternion = 1; 00561 Omega = (float)rotation["q"] ; 00562 e0 = cos(Omega*M_PI/360.0f); 00563 e1 = sin(Omega*M_PI/360.0f)* (float)rotation["n1"]; 00564 e2 = sin(Omega*M_PI/360.0f)* (float)rotation["n2"]; 00565 e3 = sin(Omega*M_PI/360.0f)* (float)rotation["n3"]; 00566 } else if ( type == "matrix" ) { 00567 is_matrix = 1; 00568 matrix[0][0] = (float)rotation["m11"]; 00569 matrix[0][1] = (float)rotation["m12"]; 00570 matrix[0][2] = (float)rotation["m13"]; 00571 matrix[1][0] = (float)rotation["m21"]; 00572 matrix[1][1] = (float)rotation["m22"]; 00573 matrix[1][2] = (float)rotation["m23"]; 00574 matrix[2][0] = (float)rotation["m31"]; 00575 matrix[2][1] = (float)rotation["m32"]; 00576 matrix[2][2] = (float)rotation["m33"]; 00577 } else { 00578 // transform_type = UNKNOWN; 00579 throw InvalidStringException(euler_type, "unknown Euler Type"); 00580 } 00581 00582 float azp = fmod(az,360.0f)*M_PI/180.0f; 00583 float altp = alt*M_PI/180.0f; 00584 float phip = fmod(phi,360.0f)*M_PI/180.0f; 00585 00586 if (!is_quaternion && !is_matrix) { 00587 matrix[0][0] = cos(phip)*cos(azp) - cos(altp)*sin(azp)*sin(phip); 00588 matrix[0][1] = cos(phip)*sin(azp) + cos(altp)*cos(azp)*sin(phip); 00589 matrix[0][2] = sin(altp)*sin(phip); 00590 matrix[1][0] = -sin(phip)*cos(azp) - cos(altp)*sin(azp)*cos(phip); 00591 matrix[1][1] = -sin(phip)*sin(azp) + cos(altp)*cos(azp)*cos(phip); 00592 matrix[1][2] = sin(altp)*cos(phip); 00593 matrix[2][0] = sin(altp)*sin(azp); 00594 matrix[2][1] = -sin(altp)*cos(azp); 00595 matrix[2][2] = cos(altp); 00596 } 00597 if (is_quaternion){ 00598 matrix[0][0] = e0 * e0 + e1 * e1 - e2 * e2 - e3 * e3; 00599 matrix[0][1] = 2.0f * (e1 * e2 + e0 * e3); 00600 matrix[0][2] = 2.0f * (e1 * e3 - e0 * e2); 00601 matrix[1][0] = 2.0f * (e2 * e1 - e0 * e3); 00602 matrix[1][1] = e0 * e0 - e1 * e1 + e2 * e2 - e3 * e3; 00603 matrix[1][2] = 2.0f * (e2 * e3 + e0 * e1); 00604 matrix[2][0] = 2.0f * (e3 * e1 + e0 * e2); 00605 matrix[2][1] = 2.0f * (e3 * e2 - e0 * e1); 00606 matrix[2][2] = e0 * e0 - e1 * e1 - e2 * e2 + e3 * e3; 00607 // keep in mind matrix[0][2] is M13 gives an e0 e2 piece, etc 00608 } 00609 00610 // Apply scale if it existed previously 00611 if (scale != 1.0f) { 00612 for(int i=0; i<3; ++i) { 00613 for(int j=0; j<3; ++j) { 00614 matrix[i][j] *= scale; 00615 } 00616 } 00617 } 00618 00619 // Apply post x mirroring if it was applied previouslys 00620 if ( x_mirror ) { 00621 for(int j=0; j<3; ++j) { 00622 matrix[0][j] *= -1.0f; 00623 } 00624 } 00625 }
void Transform::set_scale | ( | const float & | scale | ) |
Set the scale.
scale | the amount to scale by |
Definition at line 930 of file transform.cpp.
References EMAN::Util::apply_precision(), ERR_LIMIT, get_scale(), InvalidValueException, and matrix.
Referenced by EMAN::FourierReconstructor::determine_slice_agreement(), get_rotation_transform(), EMAN::BackProjectionReconstructor::insert_slice(), EMAN::FourierReconstructor::insert_slice(), set_params(), set_params_inverse(), and set_pre_trans().
00930 { 00931 if (new_scale <= 0) { 00932 throw InvalidValueException(new_scale,"The scale factor in a Transform object must be positive and non zero"); 00933 } 00934 // Transform = MTSR (Mirroring, Translation, Scaling, Rotate) 00935 // So changing the scale boils down to this.... 00936 00937 float old_scale = get_scale(); 00938 00939 float n_scale = new_scale; 00940 Util::apply_precision(n_scale,ERR_LIMIT); 00941 00942 float corrected_scale = n_scale/old_scale; 00943 if ( corrected_scale != 1.0 ) { 00944 for(int i = 0; i < 3; ++i ) { 00945 for(int j = 0; j < 3; ++j ) { 00946 matrix[i][j] *= corrected_scale; 00947 } 00948 } 00949 } 00950 }
void EMAN::Transform::set_trans | ( | const Vec2f & | v | ) | [inline] |
Set the post translation component using a Vec2f.
v | the 2D translation vector |
Definition at line 200 of file transform.h.
References set_trans(), and v.
void EMAN::Transform::set_trans | ( | const Vec3f & | v | ) | [inline] |
Set the post translation component using a Vec3f.
v | the 3D translation vector |
Definition at line 195 of file transform.h.
References set_trans(), and v.
void Transform::set_trans | ( | const float & | x, | |
const float & | y, | |||
const float & | z = 0 | |||
) |
Set the post translation component.
x | the x translation | |
y | the y translation | |
z | the z translation |
Definition at line 859 of file transform.cpp.
References get_mirror(), and matrix.
Referenced by EMAN::RefineAligner::align(), EMAN::FourierReconstructor::determine_slice_agreement(), get_pre_trans(), get_pre_trans_2d(), get_rotation_transform(), EMAN::HSym::get_sym(), EMAN::BackProjectionReconstructor::insert_slice(), EMAN::FourierReconstructor::insert_slice(), refalin3d_perturb(), set_params(), set_params_inverse(), set_pre_trans(), and set_trans().
00860 { 00861 // if ( z != 0.0 ) { 00862 // validate_and_set_type(THREED); 00863 // } 00864 bool x_mirror = get_mirror(); 00865 00866 if (x_mirror) matrix[0][3] = -x; 00867 else matrix[0][3] = x; 00868 matrix[1][3] = y; 00869 matrix[2][3] = z; 00870 }
Transform Transform::tet_3_to_2 | ( | ) | [static] |
Get the transform that moves any tetrahedron generated by eman2 so that it matches the 2-2-2 (MRC, FREALIGN) convention.
Doctor Phil says: AltAngle= acos(-1/3.0)*90/pi;
Definition at line 83 of file transform.cpp.
References t.
00083 { 00084 Transform t; 00085 Dict d; 00086 d["type"] = "eman"; 00087 d["phi"] = 45.0f; 00088 d["az"] = 0.0f; 00089 d["alt"] = 54.73561f; // 3 fold to a 2 fold 00093 t.set_rotation(d); 00094 return t; 00095 }
void Transform::to_identity | ( | ) |
Force the internal matrix to become the identity.
Definition at line 165 of file transform.cpp.
References matrix.
Referenced by Transform().
00166 { 00167 // transform_type = UNKNOWN; 00168 for(int i=0; i<3; ++i) { 00169 for(int j=0; j<4; ++j) { 00170 if(i==j) { 00171 matrix[i][j] = 1; 00172 } 00173 else { 00174 matrix[i][j] = 0; 00175 } 00176 } 00177 } 00178 }
Transform a 3D vector using the internal transformation matrix.
v | a three dimensional vector to be transformed |
Definition at line 383 of file transform.h.
References transform(), and v.
00383 { 00384 // assert_consistent_type(THREED); // Transform does the assertion 00385 return transform(v[0],v[1],v[2]); 00386 }
Vec3f EMAN::Transform::transform | ( | const float & | x, | |
const float & | y, | |||
const float & | z | |||
) | const [inline] |
Transform 3D coordinates using the internal transformation matrix.
x | the x coordinate of the transformed point | |
y | the y coordinate of the transformed point | |
z | the z coordinate of the transformed point |
Definition at line 369 of file transform.h.
References matrix.
00369 { 00370 // assert_consistent_type(THREED); 00371 Vec3f ret; 00372 ret[0] = matrix[0][0] * x + matrix[0][1] * y + matrix[0][2] * z + matrix[0][3]; 00373 ret[1] = matrix[1][0] * x + matrix[1][1] * y + matrix[1][2] * z + matrix[1][3]; 00374 ret[2] = matrix[2][0] * x + matrix[2][1] * y + matrix[2][2] * z + matrix[2][3]; 00375 return ret; 00376 }
Transform a 2D vector using the internal transformation matrix.
v | a two dimensional vector to be transformed |
Definition at line 359 of file transform.h.
References transform(), and v.
Vec2f EMAN::Transform::transform | ( | const float & | x, | |
const float & | y | |||
) | const [inline] |
Transform 2D coordinates using the internal transformation matrix.
x | the x coordinate of the transformed point | |
y | the y coordinate of the transformed point |
Definition at line 346 of file transform.h.
References matrix.
Referenced by EMAN::EMData::get_rotated_clip(), EMAN::operator *(), and transform().
00346 { 00347 // assert_valid_2d(); 00348 Vec2f ret; 00349 ret[0] = matrix[0][0]*x + matrix[0][1]*y + matrix[0][3]; 00350 ret[1] = matrix[1][0]*x + matrix[1][1]*y + matrix[1][3]; 00351 return ret; 00352 }
Transform Transform::transpose | ( | ) | const |
Get the transpose of this transformation matrix.
Definition at line 1142 of file transform.cpp.
References t.
Referenced by refalin3d_perturb(), and EMAN::PDBReader::right_transform().
void Transform::transpose_inplace | ( | ) |
Get the transpose of this transformation matrix.
Definition at line 1129 of file transform.cpp.
References matrix.
01129 { 01130 float tempij; 01131 for (int i = 0; i < 3; i++) { 01132 for (int j = 0; j < i; j++) { 01133 if (i != j) { 01134 tempij= matrix[i][j]; 01135 matrix[i][j] = matrix[j][i]; 01136 matrix[j][i] = tempij; 01137 } 01138 } 01139 } 01140 }
const float Transform::ERR_LIMIT = 0.000001f [static] |
Definition at line 87 of file transform.h.
Referenced by assert_valid_2d(), get_determinant(), get_rotation(), get_scale(), get_scale_and_mirror(), get_trans(), is_identity(), EMAN::Symmetry3D::point_in_which_asym_unit(), EMAN::Util::point_is_in_triangle_2d(), and set_scale().
float EMAN::Transform::matrix[3][4] [private] |
Definition at line 429 of file transform.h.
Referenced by assert_valid_2d(), at(), copy_matrix_into_array(), get_determinant(), get_matrix(), get_matrix3_row(), get_rotation(), get_trans(), get_trans_2d(), invert(), is_identity(), operator=(), operator[](), orthogonalize(), printme(), set(), set_matrix(), set_mirror(), set_rotation(), set_scale(), set_trans(), to_identity(), transform(), Transform(), and transpose_inplace().
vector< string > Transform::permissable_2d_not_rot [static, private] |
This map is used to validate keys in the argument maps - e.g. if the type is 2d and the angle is not "alpha" then we should throw.
Definition at line 434 of file transform.h.
Referenced by detect_problem_keys(), and init_permissable_keys().
vector< string > Transform::permissable_3d_not_rot [static, private] |
Definition at line 435 of file transform.h.
Referenced by detect_problem_keys(), and init_permissable_keys().
map< string, vector< string > > Transform::permissable_rot_keys [static, private] |
Definition at line 436 of file transform.h.
Referenced by detect_problem_keys(), and init_permissable_keys().