Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members

EMAN::Transform Class Reference
[unit test in Python]

A Transform object is a somewhat specialized object designed specifically for EMAN2/Sparx storage of alignment parameters and euler orientations. More...

#include <transform.h>

List of all members.

Public Member Functions

 Transform ()
 Default constructor Internal matrix is the identity.
 Transform (const Transform &rhs)
 Copy constructor.
Transformoperator= (const Transform &that)
 Assignment operator.
bool operator== (const Transform &rhs) const
 Equality comparision operator.
bool operator!= (const Transform &rhs) const
 Unequality comparision 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

int get_nsym (const string &sym)
 get the number of symmetries associated with the given symmetry name
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.
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

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

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.
vector< string > permissable_3d_not_rot
map< string, vector< string > > permissable_rot_keys


Detailed Description

A Transform object is a somewhat specialized object designed specifically for EMAN2/Sparx storage of alignment parameters and euler orientations.

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.

Author:
David Woolford (based on Philip Baldwin's original Transform3D code)
Date:
September 2008

Definition at line 84 of file transform.h.


Constructor & Destructor Documentation

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.

Parameters:
rhs the object to be copied

Definition at line 103 of file transform.cpp.

00104 {
00105         *this = that;
00106 }

Transform::Transform const Dict d  ) 
 

Construction using a dictionary.

Parameters:
d the dictionary containing the parameters

Definition at line 129 of file transform.cpp.

References set_params(), and to_identity().

00129                                    {
00130         to_identity();
00131         set_params(d);
00132 }

Transform::Transform const float  array[12]  ) 
 

Construction using an array of floats.

Parameters:
array the array of values that will become the internal matrix. row order (3 rows of 4)

Definition at line 135 of file transform.cpp.

References matrix.

00135                                           {
00136         memcpy(matrix,array,12*sizeof(float));
00137 }

Transform::Transform const vector< float >  array  ) 
 

Construction using a vector of size 12.

Parameters:
array the array of values that will become the internal matrix. row order (3 rows of 4)

Definition at line 139 of file transform.cpp.

References set_matrix().

00140 {
00141         set_matrix(array);
00142 }

EMAN::Transform::~Transform  )  [inline]
 

Definition at line 132 of file transform.h.

00132 { }


Member Function Documentation

void Transform::assert_valid_2d  )  const [private]
 

Definition at line 1172 of file transform.cpp.

References matrix, and UnexpectedBehaviorException.

Referenced by get_rotation(), and set_rotation().

01172                                       {
01173         int rotation_error = 0;
01174         int translation_error = 0;
01175         if (fabs(matrix[2][0]) > ERR_LIMIT) rotation_error++;
01176         if (fabs(matrix[2][1]) > ERR_LIMIT) rotation_error++;
01177         if (fabs(matrix[2][3]) > ERR_LIMIT) translation_error++;
01178         if (fabs(matrix[0][2]) > ERR_LIMIT) rotation_error++;
01179         if (fabs(matrix[1][2]) > ERR_LIMIT) rotation_error++;
01180 //      if (fabs(matrix[2][2]-1.0) >ERR_LIMIT) rotation_error++; 
01181         if ( translation_error && rotation_error ) {
01182                 throw UnexpectedBehaviorException("Error, the internal matrix contains 3D rotations and 3D translations. This object can not be considered 2D");
01183         } else if ( translation_error ) {
01184                 throw UnexpectedBehaviorException("Error, the internal matrix contains a non zero z component for a 3D translation. This object can not be considered 2D");
01185         }
01186         else if ( rotation_error ) {
01187                 throw UnexpectedBehaviorException("Error, the internal matrix contains 3D rotations and this object can not be considered 2D");
01188         }
01189 
01190 }

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 338 of file transform.h.

References matrix().

Referenced by EMAN::PawelProjector::backproject3d().

00338 { return matrix[r][c]; }

void Transform::copy_matrix_into_array float *  const  )  const
 

Definition at line 155 of file transform.cpp.

References matrix.

Referenced by EMAN::FourierReconstructor::do_compare_slice_work(), EMAN::FourierReconstructor::do_insert_slice_work(), EMAN::TransformProcessor::process(), EMAN::TransformProcessor::process_inplace(), and EMAN::StandardProjector::project3d().

00155                                                                {
00156 
00157         int idx = 0;
00158         for(int i=0; i<3; ++i) {
00159                 for(int j=0; j<4; ++j) {
00160                         array[idx] = matrix[i][j];
00161                         idx ++;
00162                 }
00163         }
00164 }

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.

Parameters:
d the dictionary that was the function argument of the set_params, set_rotation or the set_params_inv function
Exceptions:
InvalidParameterException if the dictionary is invalid in anyway

Definition at line 325 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().

00325                                                  {
00326         if (permissable_rot_keys.size() == 0 ) {
00327                 init_permissable_keys();
00328         }
00329 
00330         vector<string> verification;
00331         vector<string> problem_keys;
00332         bool is_2d = false;
00333         if (d.has_key_ci("type") ) {
00334                 string type = Util::str_to_lower((string)d["type"]);
00335                 bool problem = false;
00336                 if (permissable_rot_keys.find(type) == permissable_rot_keys.end() ) {
00337                         problem_keys.push_back(type);
00338                         problem = true;
00339                 }
00340                 if ( !problem ) {
00341                         vector<string> perm = permissable_rot_keys[type];
00342                         std::copy(perm.begin(),perm.end(),back_inserter(verification));
00343 
00344                         if ( type == "2d" ) {
00345                                 is_2d = true;
00346                                 std::copy(permissable_2d_not_rot.begin(),permissable_2d_not_rot.end(),back_inserter(verification));
00347                         }
00348                 }
00349         }
00350         if ( !is_2d ) {
00351                 std::copy(permissable_3d_not_rot.begin(),permissable_3d_not_rot.end(),back_inserter(verification));
00352         }
00353 
00354         for (Dict::const_iterator it = d.begin(); it != d.end();  ++it) {
00355                 if ( std::find(verification.begin(),verification.end(), it->first) == verification.end() ) {
00356                         problem_keys.push_back(it->first);
00357                 }
00358         }
00359 
00360         if (problem_keys.size() != 0 ) {
00361                 string error;
00362                 if (problem_keys.size() == 1) {
00363                         error = "Transform Error: The \"" +problem_keys[0]+ "\" key is unsupported";
00364                 } else {
00365                         error = "Transform Error: The ";
00366                         for(vector<string>::const_iterator cit = problem_keys.begin(); cit != problem_keys.end(); ++cit ) {
00367                                 if ( cit != problem_keys.begin() ) {
00368                                         if (cit == (problem_keys.end() -1) ) error += " and ";
00369                                         else error += ", ";
00370                                 }
00371                                 error += "\"";
00372                                 error += *cit;
00373                                 error += "\"";
00374                         }
00375                         error += " keys are unsupported";
00376                 }
00377                 throw InvalidParameterException(error);
00378         }
00379 }

float Transform::get_determinant  )  const
 

Get the determinant of the matrix.

Returns:
the determinant

Definition at line 1085 of file transform.cpp.

References EMAN::Util::apply_precision(), ERR_LIMIT, and matrix.

Referenced by get_mirror(), get_scale(), and get_scale_and_mirror().

01086 {
01087         float det;
01088         double det2;
01089         det2  = matrix[0][0]*((double)matrix[1][1]*matrix[2][2]-(double)matrix[2][1]*matrix[1][2]);
01090         det2 -= matrix[0][1]*((double)matrix[1][0]*matrix[2][2]-(double)matrix[2][0]*matrix[1][2]);
01091         det2 += matrix[0][2]*((double)matrix[1][0]*matrix[2][1]-(double)matrix[2][0]*matrix[1][1]);
01092 
01093         det = (float)det2;
01094         Util::apply_precision(det,ERR_LIMIT);
01095 
01096         return det;
01097 }

Transform Transform::get_hflip_transform  )  const
 

How do I get the transform that will yield the horizontally flipped projection?

Returns:
the transform that will yield the horizontally flipped projection

Definition at line 681 of file transform.cpp.

References get_rotation(), get_trans(), set_rotation(), set_trans(), and EMAN::Vec3f.

00681                                                {
00682 
00683         Dict rot = get_rotation("eman");
00684         rot["alt"] = 180.0f + static_cast<float>(rot["alt"]);
00685         rot["phi"] = 180.0f - static_cast<float>(rot["phi"]);
00686 
00687         Transform ret(*this); // Is the identity
00688         ret.set_rotation(rot);
00689 
00690         Vec3f trans = get_trans();
00691         trans[0] = -trans[0];
00692         ret.set_trans(trans);
00693 
00694 //      ret.set_mirror(self.get_mirror());
00695 
00696         return ret;
00697 }

vector< float > Transform::get_matrix  )  const
 

Get the transformation matrix using a vector.

Returns:
a vector - 3 rows of 4 - that stores the values of the transformation matrix

Definition at line 166 of file transform.cpp.

References matrix.

00167 {
00168         vector<float> ret(12);
00169         for(int i=0; i<3; ++i) {
00170                 for(int j=0; j<4; ++j) {
00171                         ret[i*4+j] = matrix[i][j];
00172                 }
00173         }
00174         return ret;
00175 }

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.

Parameters:
i the row number (starting at 0)
Returns:
the ith row

Definition at line 407 of file transform.h.

References matrix(), and EMAN::Vec3f.

Referenced by EMAN::PawelProjector::project3d().

00407                                                                   {
00408                                 return Vec3f(matrix[i][0], matrix[i][1], matrix[i][2]);
00409                         }

bool Transform::get_mirror  )  const
 

Query whether x_mirroring is occuring.

Returns:
whether x_mirroring is occuring

Definition at line 1056 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().

01056                                  {
01057         float determinant = get_determinant();
01058 
01059         bool x_mirror = false;
01060         if ( determinant < 0 ) x_mirror = true;
01061 
01062         return x_mirror;
01063 
01064 }

int Transform::get_nsym const string &  sym  )  [static]
 

get the number of symmetries associated with the given symmetry name

Definition at line 1203 of file transform.cpp.

References EMAN::Factory< T >::get(), and EMAN::Symmetry3D::get_nsym().

Referenced by EMAN::PointArray::set_from(), EMAN::nnSSNR_ctfReconstructor::setup(), EMAN::nn4_ctf_rectReconstructor::setup(), EMAN::nn4_ctfReconstructor::setup(), EMAN::nnSSNR_Reconstructor::setup(), EMAN::nn4_rectReconstructor::setup(), and EMAN::nn4Reconstructor::setup().

01204 {
01205         Symmetry3D* sym = Factory<Symmetry3D>::get(sym_name);
01206         int nsym = sym->get_nsym();
01207         delete sym;
01208         return nsym;
01209 }

Dict Transform::get_params const string &  euler_type  )  const
 

Get the parameters of the entire transform, using a specific euler convention.

Parameters:
euler_type the euler type of the retrieved rotation
Returns:
a dictionary containing the parameters

Definition at line 435 of file transform.cpp.

References get_mirror(), get_rotation(), get_scale(), get_trans(), EMAN::Util::str_to_lower(), v, and EMAN::Vec3f.

Referenced by EMAN::RefineAligner::align(), EMAN::Util::BPCQ(), EMAN::Util::multiref_polar_ali_2d_local(), EMAN::Util::multiref_polar_ali_2d_local_psi(), EMAN::TomoTiltEdgeMaskProcessor::process_inplace(), and EMAN::SpiderIO::write_single_header().

00435                                                          {
00436         Dict params = get_rotation(euler_type);
00437 
00438         Vec3f v = get_trans();
00439         params["tx"] = v[0]; params["ty"] = v[1];
00440 
00441         string type = Util::str_to_lower(euler_type);
00442         if ( type != "2d") params["tz"] = v[2];
00443 
00444         float scale = get_scale();
00445         params["scale"] = scale;
00446 
00447         bool mirror = get_mirror();
00448         params["mirror"] = mirror;
00449 
00450         return params;
00451 }

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.

Parameters:
euler_type the euler type of the retrieved rotation
Returns:
a dictionary containing the parameters

Definition at line 455 of file transform.cpp.

References get_mirror(), get_pre_trans(), get_rotation(), get_scale(), inverse(), EMAN::Util::str_to_lower(), v, and EMAN::Vec3f.

00455                                                                  {
00456         Transform inv(inverse());
00457 
00458         Dict params = inv.get_rotation(euler_type);
00459         Vec3f v = inv.get_pre_trans();
00460         params["tx"] = v[0]; params["ty"] = v[1];
00461 
00462         string type = Util::str_to_lower(euler_type);
00463         if ( type != "2d") params["tz"] = v[2];
00464 
00465         float scale = inv.get_scale();
00466         params["scale"] = scale;
00467 
00468         bool mirror = inv.get_mirror();
00469         params["mirror"] = mirror;
00470 
00471         return params;
00472 }

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.

Returns:
the pre translation vector

Definition at line 914 of file transform.cpp.

References get_trans(), invert(), set_trans(), and EMAN::Vec3f.

Referenced by get_params_inverse().

00915 {
00916         Transform T(*this);
00917         T.set_trans(0,0,0);
00918         T.invert();
00919 
00920         Transform soln  = T*(*this);
00921 //      soln.printme();
00922         return soln.get_trans();
00923 }

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

Returns:
the pre translation vector

Definition at line 925 of file transform.cpp.

References get_trans_2d(), invert(), set_trans(), and EMAN::Vec2f.

00926 {
00927         Transform T(*this);
00928         T.set_trans(0,0,0);
00929         T.invert();
00930 
00931         Transform soln  = T*(*this);
00932 //      soln.printme();
00933         return soln.get_trans_2d();
00934 }

Dict Transform::get_rotation const string &  euler_type = "eman"  )  const
 

Get a rotation in any Euler format.

Parameters:
euler_type the requested Euler type
Returns:
a dictionary containing the key-entry pairs describing the rotations in terms of the requested Euler type

Definition at line 715 of file transform.cpp.

References assert_valid_2d(), get_scale_and_mirror(), InvalidStringException, matrix(), matrix, phi, sqrt(), EMAN::Util::str_to_lower(), and UnexpectedBehaviorException.

Referenced by EMAN::file_store::add_image(), EMAN::RotateTranslateAligner::align(), EMAN::RotationalAligner::align(), EMAN::PawelProjector::backproject3d(), EMAN::ChaoProjector::backproject3d(), get_hflip_transform(), get_params(), get_params_inverse(), get_vflip_transform(), main(), EMAN::ChaoProjector::project3d(), EMAN::FourierGriddingProjector::project3d(), set_pre_trans(), and EMAN::RT3DSymmetryAligner::xform_align_nbest().

00716 {
00717         Dict result;
00718 
00719         //float max = 1 - ERR_LIMIT;
00720         float scale;
00721         bool x_mirror;
00722         get_scale_and_mirror(scale,x_mirror);
00723         if (scale == 0) throw UnexpectedBehaviorException("The determinant of the Transform is 0. This is unexpected.");
00724 
00725         double cosalt = matrix[2][2]/scale;
00726         double x_mirror_scale = (x_mirror ? -1.0f : 1.0f);
00727         double inv_scale = 1.0f/scale;
00728 
00729         double az  = 0;
00730         double alt = 0;
00731         double phi = 0;
00732         double phiS = 0;  // like az  (but in SPIDER ZYZ)
00733         double psiS = 0;  // like phi (but in SPIDER ZYZ)
00734 
00735         // get alt, az, phi in EMAN convention
00736 
00737         if (cosalt >= 1) {  // that is, alt close to 0
00738                         alt = 0;
00739                         az = 0;
00740                         phi = (double)EMConsts::rad2deg * atan2(x_mirror_scale*matrix[0][1], x_mirror_scale*matrix[0][0]);
00741         } else if (cosalt <= -1) {  // that is, alt close to 180
00742                         alt = 180;
00743                         az = 0;
00744                         phi = (double)EMConsts::rad2deg * atan2(-x_mirror_scale*matrix[0][1], x_mirror_scale*matrix[0][0]);
00745         } else {   // for non exceptional cases:  0 < alt < 180
00746 
00747                 az  = (double)EMConsts::rad2deg * atan2(scale*matrix[2][0], -scale*matrix[2][1]);
00748 
00749                 if (matrix[2][2]==0.0)
00750                         alt = 90.0;
00751                 else
00752                         alt = (double)EMConsts::rad2deg * atan(sqrt((double)matrix[2][0]*matrix[2][0]+(double)matrix[2][1]*matrix[2][1])/fabs(matrix[2][2]));
00753 
00754                 if (matrix[2][2] * scale < 0)
00755                         alt = 180.0f-alt;
00756                 
00757                 phi = (double)EMConsts::rad2deg * atan2(x_mirror_scale*(double)matrix[0][2], (double)matrix[1][2]);
00758 
00759         } // ends separate cases: alt close to 0, 180, or neither
00760 
00761         phi = phi-360.0*floor(phi/360.0);
00762         az  = az -360.0*floor(az/360.0);
00763 
00764 //  get phiS, psiS (SPIDER)
00765         if (cosalt >= 1) {  // that is, alt close to 0
00766                 phiS = 0;
00767                 psiS = phi;
00768         } else if (cosalt <= -1) {  // that is, alt close to 180
00769                 phiS = 0;
00770                 psiS = phi + 180.0;
00771         } else {
00772                 phiS = az  - 90.0;
00773                 psiS = phi + 90.0;
00774         }
00775 
00776         phiS = phiS-360.0*floor(phiS/360.0);
00777         psiS = psiS-360.0*floor(psiS/360.0);
00778 
00779 //   do some quaternionic stuff here
00780 
00781         double nphi = (az-phi)/2.0;
00782     // The next is also e0
00783         double cosOover2 = cos((az+phi)*EMConsts::deg2rad/2.0) * cos(alt*EMConsts::deg2rad/2.0);
00784         double sinOover2 = sqrt(1 -cosOover2*cosOover2);
00785         double cosnTheta = sin((az+phi)*EMConsts::deg2rad/2.0) * cos(alt*EMConsts::deg2rad/2.0) / sqrt(1-cosOover2*cosOover2);
00786         double sinnTheta = sqrt(1-cosnTheta*cosnTheta);
00787         double n1 = sinnTheta*cos(nphi*EMConsts::deg2rad);
00788         double n2 = sinnTheta*sin(nphi*EMConsts::deg2rad);
00789         double n3 = cosnTheta;
00790         double xtilt = 0;
00791         double ytilt = 0;
00792         double ztilt = 0;
00793 
00794 
00795         if (cosOover2<0) {
00796                 cosOover2*=-1; n1 *=-1; n2*=-1; n3*=-1;
00797         }
00798 
00799         string type = Util::str_to_lower(euler_type);
00800 
00801         result["type"] = type;
00802         if (type == "2d") {
00803                 assert_valid_2d();
00804                 result["alpha"]  = phi;
00805         } else if (type == "eman") {
00806 //              assert_consistent_type(THREED);
00807                 result["az"]  = az;
00808                 result["alt"] = alt;
00809                 result["phi"] = phi;
00810         } else if (type == "imagic") {
00811 //              assert_consistent_type(THREED);
00812                 result["alpha"] = az;
00813                 result["beta"]  = alt;
00814                 result["gamma"] = phi;
00815         } else if (type == "spider") {
00816 //              assert_consistent_type(THREED);
00817                 result["phi"]   = phiS;  // The first Euler like az
00818                 result["theta"] = alt;
00819                 result["psi"]   = psiS;
00820         } else if (type == "mrc") {
00821 //              assert_consistent_type(THREED);
00822                 result["phi"]   = phiS;
00823                 result["theta"] = alt;
00824                 result["omega"] = psiS;
00825         } else if (type == "xyz") {               // need to double-check these 3 equations ********
00826 //              assert_consistent_type(THREED);
00827                 xtilt = atan2(-sin(EMConsts::deg2rad*phiS)*sin(EMConsts::deg2rad*alt),cos(EMConsts::deg2rad*alt));
00828                 ytilt = asin(  cos(EMConsts::deg2rad*phiS)*sin(EMConsts::deg2rad*alt));
00829                 ztilt = psiS*EMConsts::deg2rad - atan2(sin(xtilt), cos(xtilt) *sin(ytilt));
00830 
00831                 xtilt *= EMConsts::rad2deg; ytilt *= EMConsts::rad2deg; ztilt *= EMConsts::rad2deg;
00832                 xtilt = xtilt-360*.0*floor((xtilt+180.0)/360.0);
00833                 ytilt = ytilt-360*.0*floor((ytilt+180.0)/360.0);  //already in range [-90,90] but anyway...
00834                 ztilt = ztilt-360*.0*floor((ztilt+180.0)/360.0);
00835 
00836                 result["xtilt"]  = xtilt;
00837                 result["ytilt"]  = ytilt;
00838                 result["ztilt"]  = ztilt;
00839         } else if (type == "quaternion") {
00840 //              assert_consistent_type(THREED);
00841                 result["e0"] = cosOover2 ;
00842                 result["e1"] = sinOover2 * n1 ;
00843                 result["e2"] = sinOover2 * n2;
00844                 result["e3"] = sinOover2 * n3;
00845         } else if (type == "spin") {
00846 //              assert_consistent_type(THREED);
00847                 result["Omega"] = 2.0*EMConsts::rad2deg * acos(cosOover2);
00848                 result["n1"] = n1;
00849                 result["n2"] = n2;
00850                 result["n3"] = n3;
00851         } else if (type == "sgirot") {
00852 //              assert_consistent_type(THREED);
00853                 result["q"] = 2.0*EMConsts::rad2deg * acos(cosOover2);
00854                 result["n1"] = n1;
00855                 result["n2"] = n2;
00856                 result["n3"] = n3;
00857         } else if (type == "matrix") {
00858 //              assert_consistent_type(THREED);
00859                 result["m11"] = x_mirror_scale*matrix[0][0]*inv_scale;
00860                 result["m12"] = x_mirror_scale*matrix[0][1]*inv_scale;
00861                 result["m13"] = x_mirror_scale*matrix[0][2]*inv_scale;
00862                 result["m21"] = matrix[1][0]*inv_scale;
00863                 result["m22"] = matrix[1][1]*inv_scale;
00864                 result["m23"] = matrix[1][2]*inv_scale;
00865                 result["m31"] = matrix[2][0]*inv_scale;
00866                 result["m32"] = matrix[2][1]*inv_scale;
00867                 result["m33"] = matrix[2][2]*inv_scale;
00868         } else {
00869                 throw InvalidStringException(euler_type, "unknown Euler Type");
00870         }
00871 
00872         return result;
00873 }

Transform Transform::get_rotation_transform  )  const
 

Get the rotation part of the tranformation matrix as a Transform object.

Returns:
a Transform describing the rotation only

Definition at line 638 of file transform.cpp.

References set_mirror(), set_scale(), and set_trans().

Referenced by EMAN::GaussFFTProjector::project3d().

00639 {
00640         Transform ret(*this);
00641         ret.set_scale(1.0);
00642         ret.set_mirror(false);
00643         ret.set_trans(0,0,0);
00644         //ret.orthogonalize(); // ?
00645         return ret;
00646 }

float Transform::get_scale  )  const
 

Get the scale that was applied.

Returns:
the scale factor

Definition at line 959 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::TransformProcessor::process(), EMAN::TransformProcessor::process_inplace(), EMAN::GaussFFTProjector::project3d(), set_pre_trans(), and set_scale().

00959                                  {
00960         float determinant = get_determinant();
00961         if (determinant < 0 ) determinant *= -1;
00962 
00963         float scale = std::pow(determinant,1.0f/3.0f);
00964         int int_scale = static_cast<int>(scale);
00965         float scale_residual = scale-static_cast<float>(int_scale);
00966         if  ( scale_residual < ERR_LIMIT ) { scale = static_cast<float>(int_scale); };
00967 
00968         Util::apply_precision(scale, ERR_LIMIT);
00969 
00970         return scale;
00971 }

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

Parameters:
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 1066 of file transform.cpp.

References EMAN::Util::apply_precision(), ERR_LIMIT, and get_determinant().

Referenced by get_rotation(), orthogonalize(), and set_rotation().

01066                                                                        {
01067 
01068         float determinant = get_determinant();
01069         x_mirror = false;
01070         if ( determinant < 0 ) {
01071                 x_mirror = true;
01072                 determinant *= -1;
01073         }
01074         if (determinant != 1 ) {
01075                 scale = std::pow(determinant,1.0f/3.0f);
01076                 int int_scale = static_cast<int>(scale);
01077                 float scale_residual = scale-static_cast<float>(int_scale);
01078                 if  ( scale_residual < ERR_LIMIT ) { scale = static_cast<float>(int_scale); };
01079         }
01080         else scale = 1;
01081 
01082         Util::apply_precision(scale,ERR_LIMIT);
01083 }

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 1194 of file transform.cpp.

References EMAN::Factory< T >::get(), and EMAN::Symmetry3D::get_sym().

Referenced by EMAN::nnSSNR_ctfReconstructor::insert_padfft_slice(), EMAN::nn4_ctfReconstructor::insert_padfft_slice(), EMAN::nnSSNR_Reconstructor::insert_padfft_slice(), EMAN::nn4_rectReconstructor::insert_padfft_slice(), EMAN::nn4Reconstructor::insert_padfft_slice(), EMAN::PointArray::set_from(), and EMAN::EMData::symvol().

01195 {
01196         Symmetry3D* sym = Factory<Symmetry3D>::get(sym_name);
01197         Transform ret;
01198         ret = (*this) * sym->get_sym(n);
01199         delete sym;
01200         return ret;
01201 }

Vec3f Transform::get_trans  )  const
 

Get the post trans as a vec3f.

Returns:
the translation vector

Definition at line 885 of file transform.cpp.

References EMAN::Util::apply_precision(), ERR_LIMIT, get_mirror(), matrix, v, and EMAN::Vec3f.

Referenced by EMAN::Refine3DAlignerGrid::align(), get_hflip_transform(), get_params(), get_pre_trans(), get_vflip_transform(), EMAN::ACFCenterProcessor::process_inplace(), EMAN::GaussFFTProjector::project3d(), EMAN::EMData::rot_scale_trans(), EMAN::EMData::rot_scale_trans_background(), and set_params_inverse().

00886 {
00887         // No type asserted
00888         bool x_mirror = get_mirror();
00889         Vec3f v;
00890         if (x_mirror) v[0] = -matrix[0][3];
00891         else v[0] = matrix[0][3];
00892         v[1] = matrix[1][3];
00893         v[2] = matrix[2][3];
00894 
00895         Util::apply_precision(v[0],ERR_LIMIT);
00896         Util::apply_precision(v[1],ERR_LIMIT);
00897         Util::apply_precision(v[2],ERR_LIMIT);
00898 
00899         return v;
00900 }

Vec2f Transform::get_trans_2d  )  const
 

Get the degenerant 2D post trans as a vec2f.

Returns:
the 2D translation vector

Definition at line 902 of file transform.cpp.

References get_mirror(), matrix, v, and EMAN::Vec2f.

Referenced by get_pre_trans_2d(), EMAN::padfft_slice(), and EMAN::BackProjectionReconstructor::preprocess_slice().

00903 {
00904         bool x_mirror = get_mirror();
00905         Vec2f v;
00906         if (x_mirror) v[0] = -matrix[0][3];
00907         else v[0] = matrix[0][3];
00908         v[1] = matrix[1][3];
00909         return v;
00910 }

Transform Transform::get_vflip_transform  )  const
 

How do I get the transform that will yield the vertically flipped projection?

Returns:
the transform that will yield the vertically flipped projection

Definition at line 699 of file transform.cpp.

References get_rotation(), get_trans(), set_rotation(), set_trans(), and EMAN::Vec3f.

00699                                                {
00700 
00701         Dict rot = get_rotation("eman");
00702         rot["alt"] = 180.0f + static_cast<float>(rot["alt"]);
00703         rot["phi"] = - static_cast<float>(rot["phi"]);
00704 
00705         Transform ret(*this);
00706         ret.set_rotation(rot);
00707 
00708         Vec3f trans = get_trans();
00709         trans[1] = -trans[1];
00710         ret.set_trans(trans);
00711 
00712         return ret;
00713 }

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.

Returns:
a Transforms, alt=58.282523, az=270

Definition at line 67 of file transform.cpp.

References set_rotation(), and 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 240 of file transform.cpp.

References permissable_2d_not_rot, permissable_3d_not_rot, and permissable_rot_keys.

Referenced by detect_problem_keys().

00241 {
00242 
00243         permissable_2d_not_rot.push_back("tx");
00244         permissable_2d_not_rot.push_back("ty");
00245         permissable_2d_not_rot.push_back("scale");
00246         permissable_2d_not_rot.push_back("mirror");
00247         permissable_2d_not_rot.push_back("type");
00248 
00249         permissable_3d_not_rot.push_back("tx");
00250         permissable_3d_not_rot.push_back("ty");
00251         permissable_3d_not_rot.push_back("tz");
00252         permissable_3d_not_rot.push_back("scale");
00253         permissable_3d_not_rot.push_back("mirror");
00254         permissable_3d_not_rot.push_back("type");
00255 
00256         vector<string> tmp;
00257         tmp.push_back("alpha");
00258         permissable_rot_keys["2d"] = tmp;
00259 
00260         tmp.clear();
00261         tmp.push_back("alt");
00262         tmp.push_back("az");
00263         tmp.push_back("phi");
00264         permissable_rot_keys["eman"] = tmp;
00265 
00266         tmp.clear();
00267         tmp.push_back("psi");
00268         tmp.push_back("theta");
00269         tmp.push_back("phi");
00270         permissable_rot_keys["spider"] = tmp;
00271 
00272         tmp.clear();
00273         tmp.push_back("alpha");
00274         tmp.push_back("beta");
00275         tmp.push_back("gamma");
00276         permissable_rot_keys["imagic"] = tmp;
00277 
00278         tmp.clear();
00279         tmp.push_back("ztilt");
00280         tmp.push_back("xtilt");
00281         tmp.push_back("ytilt");
00282         permissable_rot_keys["xyz"] = tmp;
00283 
00284         tmp.clear();
00285         tmp.push_back("phi");
00286         tmp.push_back("theta");
00287         tmp.push_back("omega");
00288         permissable_rot_keys["mrc"] = tmp;
00289 
00290         tmp.clear();
00291         tmp.push_back("e0");
00292         tmp.push_back("e1");
00293         tmp.push_back("e2");
00294         tmp.push_back("e3");
00295         permissable_rot_keys["quaternion"] = tmp;
00296 
00297         tmp.clear();
00298         tmp.push_back("n1");
00299         tmp.push_back("n2");
00300         tmp.push_back("n3");
00301         tmp.push_back("Omega");
00302         permissable_rot_keys["spin"] = tmp;
00303 
00304         tmp.clear();
00305         tmp.push_back("n1");
00306         tmp.push_back("n2");
00307         tmp.push_back("n3");
00308         tmp.push_back("q");
00309         permissable_rot_keys["sgirot"] = tmp;
00310 
00311         tmp.clear();
00312         tmp.push_back("m11");
00313         tmp.push_back("m12");
00314         tmp.push_back("m13");
00315         tmp.push_back("m21");
00316         tmp.push_back("m22");
00317         tmp.push_back("m23");
00318         tmp.push_back("m31");
00319         tmp.push_back("m32");
00320         tmp.push_back("m33");
00321         permissable_rot_keys["matrix"] = tmp;
00322 }

Transform Transform::inverse  )  const
 

Get the inverse of this transformation matrix.

Returns:
the inverse of this transformation matrix

Definition at line 1133 of file transform.cpp.

References invert(), and t.

Referenced by get_params_inverse(), EMAN::EMData::max_3D_pixel_error(), EMAN::TransformProcessor::process(), EMAN::TransformProcessor::process_inplace(), EMAN::StandardProjector::project3d(), EMAN::EMData::rot_scale_trans(), EMAN::EMData::rot_scale_trans_background(), and EMAN::TransformProcessor::transform().

01133                                    {
01134         Transform t(*this);
01135         t.invert();
01136         return t;
01137 }

void Transform::invert  ) 
 

Get the inverse of this transformation matrix.

Returns:
the inverse of this transformation matrix

Definition at line 1099 of file transform.cpp.

References matrix, and v0.

Referenced by EMAN::EMData::extract_plane(), EMAN::EMData::extract_plane_rect(), EMAN::EMData::extract_plane_rect_fast(), get_pre_trans(), get_pre_trans_2d(), EMAN::Symmetry3D::in_which_asym_unit(), inverse(), EMAN::GaussFFTProjector::project3d(), EMAN::Symmetry3D::reduce(), set_params_inverse(), set_pre_trans(), and EMAN::RT3DSphereAligner::xform_align_nbest().

01099                        {
01100 
01101         double m00 = matrix[0][0]; double m01=matrix[0][1]; double m02=matrix[0][2];
01102         double m10 = matrix[1][0]; double m11=matrix[1][1]; double m12=matrix[1][2];
01103         double m20 = matrix[2][0]; double m21=matrix[2][1]; double m22=matrix[2][2];
01104         double v0  = matrix[0][3]; double v1 =matrix[1][3]; double v2 =matrix[2][3];
01105 
01106         double cof00 = m11*m22-m12*m21;
01107         double cof11 = m22*m00-m20*m02;
01108         double cof22 = m00*m11-m01*m10;
01109         double cof01 = m10*m22-m20*m12;
01110         double cof02 = m10*m21-m20*m11;
01111         double cof12 = m00*m21-m01*m20;
01112         double cof10 = m01*m22-m02*m21;
01113         double cof20 = m01*m12-m02*m11;
01114         double cof21 = m00*m12-m10*m02;
01115 
01116         double det = m00* cof00 + m02* cof02 -m01*cof01;
01117 
01118         matrix[0][0] =   (float)(cof00/det);
01119         matrix[0][1] = - (float)(cof10/det);
01120         matrix[0][2] =   (float)(cof20/det);
01121         matrix[1][0] = - (float)(cof01/det);
01122         matrix[1][1] =   (float)(cof11/det);
01123         matrix[1][2] = - (float)(cof21/det);
01124         matrix[2][0] =   (float)(cof02/det);
01125         matrix[2][1] = - (float)(cof12/det);
01126         matrix[2][2] =   (float)(cof22/det);
01127 
01128         matrix[0][3] =  (float)((- cof00*v0 + cof10*v1 - cof20*v2)/det);
01129         matrix[1][3] =  (float)((  cof01*v0 - cof11*v1 + cof21*v2)/det);
01130         matrix[2][3] =  (float)((- cof02*v0 + cof12*v1 - cof22*v2)/det);
01131 }

bool Transform::is_identity  )  const
 

Returns whethers or this matrix is the identity.

Definition at line 192 of file transform.cpp.

References EMAN::Util::apply_precision(), ERR_LIMIT, and matrix.

Referenced by EMAN::TestTomoImage::insert_rectangle(), and EMAN::FourierReconstructor::preprocess_slice().

00192                                   {
00193         for(int i=0; i<3; ++i) {
00194                 for(int j=0; j<4; ++j) {
00195                         float c = matrix[i][j];
00196                         Util::apply_precision(c,ERR_LIMIT);
00197                         if(i==j) {
00198                                 if (c != 1.0) return false;
00199                         }
00200                         else {
00201                                 if (c != 0.0) return false;
00202                         }
00203                 }
00204         }
00205         return true;
00206 }

Transform Transform::negate  )  const
 

Negates the Transform - a useful way, for example, for getting an orientation on the opposite side of the sphere.

Returns:
a transform that has been negated

Definition at line 670 of file transform.cpp.

References set(), and t.

00671 {
00672         Transform t(*this);
00673         for(unsigned int i = 0; i < 3; ++i) {
00674                 for(unsigned int j = 0; j < 4; ++j) {
00675                         t.set(i,j,t[i][j]*-1);
00676                 }
00677         }
00678         return t;
00679 }

bool Transform::operator!= const Transform rhs  )  const
 

Unequality comparision operator.

Parameters:
rhs the Transform object compared to
Returns:
a boolean indicate if the pass in transform is not equal this object

Definition at line 125 of file transform.cpp.

References operator==(), and rhs.

00125                                                     {
00126         return !(operator==(rhs));
00127 }

Transform & Transform::operator= const Transform that  ) 
 

Assignment operator.

Parameters:
that that which this will become

Definition at line 108 of file transform.cpp.

References matrix.

00108                                                       {
00109         if (this != &that ) {
00110                 memcpy(matrix,that.matrix,12*sizeof(float));
00111 //              transform_type = that.transform_type;
00112         }
00113         return *this;
00114 }

bool Transform::operator== const Transform rhs  )  const
 

Equality comparision operator.

Parameters:
rhs the Transform object compared to
Returns:
a boolean indicate if the pass in transform is equal this object

Definition at line 116 of file transform.cpp.

References matrix, and rhs.

Referenced by operator!=().

00116                                                     {
00117         if (memcmp(this->matrix, rhs.matrix, 3*4*sizeof(float)) == 0) {
00118                 return true;
00119         }
00120         else {
00121                 return false;
00122         }
00123 }

const float* EMAN::Transform::operator[] int  i  )  const [inline]
 

Operator[] convenience so Transform3D[2][2] etc terminology can be used.

Definition at line 352 of file transform.h.

References matrix().

00352 { return matrix[i]; }

float* EMAN::Transform::operator[] int  i  )  [inline]
 

Operator[] convenience so Transform3D[2][2] etc terminology can be used.

Definition at line 347 of file transform.h.

References matrix().

00347 { 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 985 of file transform.cpp.

References get_scale_and_mirror(), matrix, R, UnexpectedBehaviorException, and V.

00986 {
00987         float scale;
00988         bool x_mirror;
00989         get_scale_and_mirror(scale,x_mirror);
00990         if (scale == 0) throw UnexpectedBehaviorException("The determinant of the Transform is 0. This is unexpected.");
00991         double inv_scale = 1.0/static_cast<double>(scale);
00992         double mirror_scale = (x_mirror == true ? -1.0:1.0);
00993 
00994         gsl_matrix * R = gsl_matrix_calloc(3,3);
00995         for ( unsigned int i = 0; i < 3; ++i )
00996         {
00997                 for ( unsigned int j = 0; j < 3; ++j )
00998                 {
00999                         if (i == 0 && mirror_scale != 1.0 ) {
01000                                 gsl_matrix_set( R, i, j, static_cast<double>(matrix[i][j])*mirror_scale*inv_scale );
01001                         }
01002                         else {
01003                                 gsl_matrix_set( R, i, j, static_cast<double>(matrix[i][j])*inv_scale );
01004                         }
01005                 }
01006         }
01007 
01008         gsl_matrix * V = gsl_matrix_calloc(3,3);
01009         gsl_vector * S = gsl_vector_calloc(3);
01010         gsl_vector * work = gsl_vector_calloc(3);
01011         gsl_linalg_SV_decomp (R, V, S, work); // Now R is U of the SVD R = USV^T
01012 
01013         gsl_matrix * Soln = gsl_matrix_calloc(3,3);
01014         gsl_blas_dgemm (CblasNoTrans, CblasTrans, 1.0, R, V, 0.0, Soln);
01015 
01016         for ( unsigned int i = 0; i < 3; ++i )
01017         {
01018                 for ( unsigned int j = 0; j < 3; ++j )
01019                 {
01020                         matrix[i][j] = static_cast<float>( gsl_matrix_get(Soln,i,j) );
01021                 }
01022         }
01023 
01024         // Apply scale if it existed previously
01025         if (scale != 1.0f) {
01026                 for(int i=0; i<3; ++i) {
01027                         for(int j=0; j<3; ++j) {
01028                                 matrix[i][j] *= scale;
01029                         }
01030                 }
01031         }
01032 
01033         // Apply post x mirroring if it was applied previouslys
01034         if ( x_mirror ) {
01035                 for(int j=0; j<3; ++j) {
01036                         matrix[0][j] *= -1.0f;
01037                 }
01038         }
01039 
01040         gsl_matrix_free(V); gsl_matrix_free(R); gsl_matrix_free(Soln);
01041         gsl_vector_free(S); gsl_vector_free(work);
01042 }

void EMAN::Transform::printme  )  const [inline]
 

Print the contents of the internal matrix verbatim to standard out.

Definition at line 298 of file transform.h.

References matrix().

00298                                              {
00299                                 printf("%8.6f %8.6f %8.6f %8.6f\n",matrix[0][0],matrix[0][1],matrix[0][2],matrix[0][3]);
00300                                 printf("%8.6f %8.6f %8.6f %8.6f\n",matrix[1][0],matrix[1][1],matrix[1][2],matrix[1][3]);
00301                                 printf("%8.6f %8.6f %8.6f %8.6f\n",matrix[2][0],matrix[2][1],matrix[2][2],matrix[2][3]);
00302                                 printf("%8.6f %8.6f %8.6f %8.6f\n",0.0,0.0,0.0,1.0);
00303 
00304                         }

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 342 of file transform.h.

References matrix().

Referenced by EMAN::PointArray::align_2d(), and negate().

00342 { matrix[r][c] = value; }

void Transform::set_matrix const vector< float > &  v  ) 
 

Set the transformation matrix using a vector.

Must be of length 12.

Parameters:
v the transformation matrix stored as a vector - 3 rows of 4.

Definition at line 144 of file transform.cpp.

References InvalidParameterException, matrix, and v.

Referenced by EMAN::EMObject::operator Transform *(), and Transform().

00145 {
00146         if (v.size() != 12 ) throw InvalidParameterException("The construction array must be of size 12");
00147 
00148         for(int i=0; i<3; ++i) {
00149                 for(int j=0; j<4; ++j) {
00150                         matrix[i][j] = v[i*4+j];
00151                 }
00152         }
00153 }

void Transform::set_mirror const bool  x_mirror  ) 
 

Set whether or not x_mirroring is occuring.

Parameters:
x_mirror whether x_mirroring should be applied

Definition at line 1044 of file transform.cpp.

References get_mirror(), and matrix.

Referenced by EMAN::RTFSlowExhaustiveAligner::align(), EMAN::RTFExhaustiveAligner::align(), EMAN::RotateFlipAlignerIterative::align(), EMAN::RotateFlipAligner::align(), EMAN::RotateTranslateFlipAlignerIterative::align(), EMAN::RotateTranslateFlipAligner::align(), EMAN::WienerFourierReconstructor::determine_slice_agreement(), EMAN::FourierReconstructor::determine_slice_agreement(), get_rotation_transform(), EMAN::WienerFourierReconstructor::insert_slice(), EMAN::FourierReconstructor::insert_slice(), refalifn(), set_params(), and set_params_inverse().

01044                                                {
01045 
01046         bool old_x_mirror = get_mirror();
01047         if (old_x_mirror == x_mirror) return; // The user is setting the same value
01048         else {
01049                 // Toggle the mirroring operation
01050                 for (int j = 0; j < 4; ++j ) {
01051                         matrix[0][j] *= -1;
01052                 }
01053         }
01054 }

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"

Parameters:
d the dictionary containing the parameters

Definition at line 209 of file transform.cpp.

References detect_problem_keys(), EMAN::Dict::get_ci(), EMAN::EMObject::get_type(), EMAN::Dict::has_key_ci(), InvalidParameterException, set_mirror(), set_rotation(), set_scale(), and set_trans().

Referenced by Transform().

00209                                         {
00210         detect_problem_keys(d);
00211 
00212         if (d.has_key_ci("type") ) set_rotation(d);
00213 
00214         if (d.has_key_ci("scale")) {
00215                 float scale = static_cast<float>(d.get_ci("scale"));
00216                 set_scale(scale);
00217         }
00218 
00219         float dx=0,dy=0,dz=0;
00220 
00221         if (d.has_key_ci("tx")) dx = static_cast<float>(d.get_ci("tx"));
00222         if (d.has_key_ci("ty")) dy = static_cast<float>(d.get_ci("ty"));
00223         if (d.has_key_ci("tz")) dz = static_cast<float>(d.get_ci("tz"));
00224 
00225         if ( dx != 0.0 || dy != 0.0 || dz != 0.0 ) {
00226                 set_trans(dx,dy,dz);
00227         }
00228 
00229         if (d.has_key_ci("mirror")) {
00230                 EMObject e = d.get_ci("mirror");
00231                 if ( (e.get_type() != EMObject::BOOL ) && (e.get_type() != EMObject::INT ) && (e.get_type() != EMObject::UNSIGNEDINT ) )
00232                         throw InvalidParameterException("Error, mirror must be a bool or an int");
00233 
00234                 bool mirror = static_cast<bool>(e);
00235                 set_mirror(mirror);
00236         }
00237 }

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.

Parameters:
d the dictionary containing the inverse parameters

Definition at line 381 of file transform.cpp.

References detect_problem_keys(), EMAN::Dict::get_ci(), get_trans(), EMAN::EMObject::get_type(), EMAN::Dict::has_key_ci(), InvalidParameterException, invert(), set_mirror(), set_rotation(), set_scale(), set_trans(), and Transform().

00381                                                 {
00382         detect_problem_keys(d);
00383 
00384         if (d.has_key_ci("type") ) set_rotation(d);
00385 
00386         float dx=0,dy=0,dz=0;
00387         if (d.has_key_ci("tx")) dx = static_cast<float>(d.get_ci("tx"));
00388         if (d.has_key_ci("ty")) dy = static_cast<float>(d.get_ci("ty"));
00389         if (d.has_key_ci("tz")) dz = static_cast<float>(d.get_ci("tz"));
00390 
00391         if ( (dx != 0.0 || dy != 0.0 || dz != 0.0) && d.has_key_ci("type") ) {
00392                 Transform pre_trans;
00393                 pre_trans.set_trans(dx,dy,dz);
00394 
00395                 Transform tmp;
00396                 tmp.set_rotation(d);
00397 
00398                 if (d.has_key_ci("scale")) {
00399                         float scale = static_cast<float>(d.get_ci("scale"));
00400                         tmp.set_scale(scale);
00401                 }
00402 
00403                 Transform solution_trans = tmp*pre_trans;
00404 
00405                 if (d.has_key_ci("scale")) {
00406                         Transform tmp;
00407                         float scale = static_cast<float>(d.get_ci("scale"));
00408                         tmp.set_scale(scale);
00409                         solution_trans = solution_trans*tmp;
00410                 }
00411 
00412                 tmp = Transform();
00413                 tmp.set_rotation(d);
00414                 solution_trans = solution_trans*tmp;
00415                 set_trans(solution_trans.get_trans());
00416         }
00417 
00418         if (d.has_key_ci("scale")) {
00419                 float scale = static_cast<float>(d.get_ci("scale"));
00420                 set_scale(scale);
00421         }
00422 
00423         if (d.has_key_ci("mirror")) {
00424                 EMObject e = d.get_ci("mirror");
00425                 if ( (e.get_type() != EMObject::BOOL ) && (e.get_type() != EMObject::INT ) && (e.get_type() != EMObject::UNSIGNEDINT ) )
00426                         throw InvalidParameterException("Error, mirror must be a bool or an int");
00427 
00428                 bool mirror = static_cast<bool>(e);
00429                 set_mirror(mirror);
00430         }
00431         invert();
00432 }

template<typename type>
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.

Parameters:
v the vector (Vec3f or Vec2f) that is the pre trans

Definition at line 495 of file transform.h.

References get_rotation(), get_scale(), invert(), set_rotation(), set_scale(), set_trans(), and v.

Referenced by EMAN::RTFExhaustiveAligner::align(), EMAN::PointArray::align_2d(), and EMAN::EMData::rotate_translate().

00495                                                    {
00496 
00497                 Transform tmp;
00498                 Dict rot = get_rotation("eman");
00499                 tmp.set_rotation(rot);
00500 
00501                 float scale = get_scale();
00502                 if (scale != 1.0 ) tmp.set_scale(scale);
00503 
00504                 Transform trans;
00505                 trans.set_trans(v);
00506 
00507                 trans = tmp*trans;
00508 
00509                 Transform tmp2;
00510                 tmp2.set_rotation(rot);
00511                 tmp2.invert(); // invert
00512                 if (scale != 1.0 ) tmp2.set_scale(1.0f/scale);
00513 
00514 
00515                 trans = trans*tmp2;
00516 
00517                 set_trans(trans.get_trans());
00518         }

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.

Parameters:
v the direction you want to solve for

Definition at line 648 of file transform.cpp.

References EMAN::Vec3< Type >::normalize(), set_rotation(), theta, UnexpectedBehaviorException, v, and EMAN::Vec3f.

00649 {
00650         if ( v[0] == 0 && v[1] == 0 && v[2] == 0 )
00651                 throw UnexpectedBehaviorException("Can't set rotation for the null vector");
00652 
00653         Vec3f v1(v);
00654         v1.normalize();
00655 
00656         double theta = acos(v1[2]); // in radians
00657         double psi = atan2(v1[1],-v1[0]);
00658 
00659         Dict d;
00660         d["theta"] = (double)EMConsts::rad2deg*theta;
00661         d["psi"] = (double)EMConsts::rad2deg*psi;
00662         d["phi"] = (double)0.0;
00663         d["type"] = "spider";
00664 
00665         set_rotation(d);
00666 
00667 
00668 }

void Transform::set_rotation const Dict rotation  ) 
 

Set a rotation using a specific Euler type and the dictionary interface Works for all Euler types.

Parameters:
rotation a dictionary containing all key-entry pair required of the associated Euler type

Definition at line 475 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::RotateTranslateAligner::align(), EMAN::RotatePrecenterAligner::align(), get_hflip_transform(), get_vflip_transform(), icos_5_to_2(), EMAN::BackProjectionReconstructor::preprocess_slice(), EMAN::FourierReconstructor::preprocess_slice(), EMAN::SymAlignProcessor::process(), EMAN::TestTomoImage::process_inplace(), EMAN::EMData::rotate_translate(), set_params(), set_params_inverse(), set_pre_trans(), set_rotation(), tet_3_to_2(), and EMAN::RT3DSphereAligner::xform_align_nbest().

00476 {
00477         detect_problem_keys(rotation);
00478         string euler_type;
00479 
00480         if (!rotation.has_key_ci("type") ){
00481                         throw InvalidParameterException("argument dictionary does not contain the type key");
00482         }
00483 
00484         euler_type = static_cast<string>(rotation.get_ci("type"));// Warning, will throw
00485 
00486 
00487         double e0=0; double e1=0; double e2=0; double e3=0;
00488         double Omega=0;
00489         double az  = 0;
00490         double alt = 0;
00491         double phi = 0;
00492         double cxtilt = 0;
00493         double sxtilt = 0;
00494         double cytilt = 0;
00495         double sytilt = 0;
00496         double cztilt = 0;
00497         double sztilt = 0;
00498         bool is_quaternion = 0;
00499         bool is_matrix = 0;
00500         bool is_xyz = 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 = (double)rotation["alpha"] ;
00514         } else if ( type == "eman" ) {
00515 //              validate_and_set_type(THREED);
00516                 az  = (double)rotation["az"] ;
00517                 alt = (double)rotation["alt"]  ;
00518                 phi = (double)rotation["phi"] ;
00519         } else if ( type == "imagic" ) {
00520 //              validate_and_set_type(THREED);
00521                 az  = (double)rotation["alpha"] ;
00522                 alt = (double)rotation["beta"]  ;
00523                 phi = (double)rotation["gamma"] ;
00524         } else if ( type == "spider" ) {
00525 //              validate_and_set_type(THREED);
00526                 az =  (double)rotation["phi"]    + 90.0;
00527                 alt = (double)rotation["theta"] ;
00528                 phi = (double)rotation["psi"]    - 90.0;
00529         } else if ( type == "xyz" ) {
00530 //              validate_and_set_type(THREED);
00531                 is_xyz = 1;
00532                 cxtilt = cos(EMConsts::deg2rad*(double)rotation["xtilt"]);
00533                 sxtilt = sin(EMConsts::deg2rad*(double)rotation["xtilt"]);
00534                 cytilt = cos(EMConsts::deg2rad*(double)rotation["ytilt"]);
00535                 sytilt = sin(EMConsts::deg2rad*(double)rotation["ytilt"]);
00536                 cztilt = cos(EMConsts::deg2rad*(double)rotation["ztilt"]);
00537                 sztilt = sin(EMConsts::deg2rad*(double)rotation["ztilt"]);
00538         } else if ( type == "mrc" ) {
00539 //              validate_and_set_type(THREED);
00540                 az  = (double)rotation["phi"]   + 90.0f ;
00541                 alt = (double)rotation["theta"] ;
00542                 phi = (double)rotation["omega"] - 90.0f ;
00543         } else if ( type == "quaternion" ) {
00544 //              validate_and_set_type(THREED);
00545                 is_quaternion = 1;
00546                 e0 = (double)rotation["e0"];
00547                 e1 = (double)rotation["e1"];
00548                 e2 = (double)rotation["e2"];
00549                 e3 = (double)rotation["e3"];
00550         } else if ( type == "spin" ) {
00551 //              validate_and_set_type(THREED);
00552                 is_quaternion = 1;
00553                 Omega = (double)rotation["Omega"];
00554                 e0 = cos(Omega*EMConsts::deg2rad/2.0);
00555                 e1 = sin(Omega*EMConsts::deg2rad/2.0) * (double)rotation["n1"];
00556                 e2 = sin(Omega*EMConsts::deg2rad/2.0) * (double)rotation["n2"];
00557                 e3 = sin(Omega*EMConsts::deg2rad/2.0) * (double)rotation["n3"];
00558         } else if ( type == "sgirot" ) {
00559 //              validate_and_set_type(THREED);
00560                 is_quaternion = 1;
00561                 Omega = (double)rotation["q"] ;
00562                 e0 = cos(Omega*EMConsts::deg2rad/2.0);
00563                 e1 = sin(Omega*EMConsts::deg2rad/2.0) * (double)rotation["n1"];
00564                 e2 = sin(Omega*EMConsts::deg2rad/2.0) * (double)rotation["n2"];
00565                 e3 = sin(Omega*EMConsts::deg2rad/2.0) * (double)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         double azp  =  az*EMConsts::deg2rad;
00583         double altp = alt*EMConsts::deg2rad;
00584         double phip = phi*EMConsts::deg2rad;
00585 
00586         if (!is_quaternion && !is_matrix && !is_xyz) {
00587                 matrix[0][0] =  (float)(cos(phip)*cos(azp) - cos(altp)*sin(azp)*sin(phip));
00588                 matrix[0][1] =  (float)(cos(phip)*sin(azp) + cos(altp)*cos(azp)*sin(phip));
00589                 matrix[0][2] =  (float)(sin(altp)*sin(phip));
00590                 matrix[1][0] =  (float)(-sin(phip)*cos(azp) - cos(altp)*sin(azp)*cos(phip));
00591                 matrix[1][1] =  (float)(-sin(phip)*sin(azp) + cos(altp)*cos(azp)*cos(phip));
00592                 matrix[1][2] =  (float)(sin(altp)*cos(phip));
00593                 matrix[2][0] =  (float)(sin(altp)*sin(azp));
00594                 matrix[2][1] =  (float)(-sin(altp)*cos(azp));
00595                 matrix[2][2] =  (float)cos(altp);
00596         }
00597         if (is_quaternion){
00598                 matrix[0][0] = (float)(e0 * e0 + e1 * e1 - e2 * e2 - e3 * e3);
00599                 matrix[0][1] = (float)(2.0f * (e1 * e2 + e0 * e3));
00600                 matrix[0][2] = (float)(2.0f * (e1 * e3 - e0 * e2));
00601                 matrix[1][0] = (float)(2.0f * (e2 * e1 - e0 * e3));
00602                 matrix[1][1] = (float)(e0 * e0 - e1 * e1 + e2 * e2 - e3 * e3);
00603                 matrix[1][2] = (float)(2.0f * (e2 * e3 + e0 * e1));
00604                 matrix[2][0] = (float)(2.0f * (e3 * e1 + e0 * e2));
00605                 matrix[2][1] = (float)(2.0f * (e3 * e2 - e0 * e1));
00606                 matrix[2][2] = (float)(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         if (is_xyz){
00610                 matrix[0][0] =  (float)(cytilt*cztilt);
00611                 matrix[0][1] =  (float)(cxtilt*sztilt+sxtilt*sytilt*cztilt);
00612                 matrix[0][2] =  (float)(sxtilt*sztilt-cxtilt*sytilt*cztilt);
00613                 matrix[1][0] =  (float)(-cytilt*sztilt);
00614                 matrix[1][1] =  (float)(cxtilt*cztilt-sxtilt*sytilt*sztilt);
00615                 matrix[1][2] =  (float)(sxtilt*cztilt+cxtilt*sytilt*sztilt);
00616                 matrix[2][0] =  (float)(sytilt);
00617                 matrix[2][1] =  (float)(-sxtilt*cytilt);
00618                 matrix[2][2] =  (float)(cxtilt*cytilt);
00619         }
00620         
00621         // Apply scale if it existed previously
00622         if (scale != 1.0f) {
00623                 for(int i=0; i<3; ++i) {
00624                         for(int j=0; j<3; ++j) {
00625                                 matrix[i][j] *= scale;
00626                         }
00627                 }
00628         }
00629 
00630         // Apply post x mirroring if it was applied previously
00631         if ( x_mirror ) {
00632                 for(int j=0; j<3; ++j) {
00633                         matrix[0][j] *= -1.0f;
00634                 }
00635         }
00636 }

void Transform::set_scale const float &  scale  ) 
 

Set the scale.

Parameters:
scale the amount to scale by

Definition at line 937 of file transform.cpp.

References EMAN::Util::apply_precision(), ERR_LIMIT, get_scale(), InvalidValueException, and matrix.

Referenced by EMAN::WienerFourierReconstructor::determine_slice_agreement(), EMAN::FourierReconstructor::determine_slice_agreement(), get_rotation_transform(), EMAN::WienerFourierReconstructor::insert_slice(), EMAN::FourierReconstructor::insert_slice(), main(), EMAN::ScaleTransformProcessor::process(), EMAN::ScaleTransformProcessor::process_inplace(), refalifn(), EMAN::EMData::scale(), set_params(), set_params_inverse(), and set_pre_trans().

00937                                                 {
00938         if (new_scale <= 0) {
00939                 throw InvalidValueException(new_scale,"The scale factor in a Transform object must be positive and non zero");
00940         }
00941         // Transform = MTSR (Mirroring, Translation, Scaling, Rotate)
00942         // So changing the scale boils down to this....
00943 
00944         float old_scale = get_scale();
00945 
00946         float n_scale = new_scale;
00947         Util::apply_precision(n_scale,ERR_LIMIT);
00948 
00949         float corrected_scale = n_scale/old_scale;
00950         if ( corrected_scale != 1.0 ) {
00951                 for(int i = 0; i < 3;  ++i ) {
00952                         for(int j = 0; j < 3; ++j ) {
00953                                 matrix[i][j] *= corrected_scale;
00954                         }
00955                 }
00956         }
00957 }

void EMAN::Transform::set_trans const Vec2f v  )  [inline]
 

Set the post translation component using a Vec2f.

Parameters:
v the 2D translation vector

Definition at line 212 of file transform.h.

References v, and EMAN::Vec2f.

00212 { set_trans(v[0],v[1]); }

void EMAN::Transform::set_trans const Vec3f v  )  [inline]
 

Set the post translation component using a Vec3f.

Parameters:
v the 3D translation vector

Definition at line 207 of file transform.h.

References v, and EMAN::Vec3f.

00207 { set_trans(v[0],v[1],v[2]); }

void Transform::set_trans const float &  x,
const float &  y,
const float &  z = 0
 

Set the post translation component.

Parameters:
x the x translation
y the y translation
z the z translation

Definition at line 875 of file transform.cpp.

References get_mirror(), and matrix.

Referenced by EMAN::Refine3DAlignerGrid::align(), EMAN::RTFSlowExhaustiveAligner::align(), EMAN::TranslationalAligner::align(), EMAN::WienerFourierReconstructor::determine_slice_agreement(), EMAN::FourierReconstructor::determine_slice_agreement(), EMAN::TestUtil::emobject_transformarray_to_py(), frm_2d_Align(), EMAN::TestUtil::get_debug_transform(), get_hflip_transform(), get_pre_trans(), get_pre_trans_2d(), get_rotation_transform(), EMAN::HSym::get_sym(), get_vflip_transform(), EMAN::WienerFourierReconstructor::insert_slice(), EMAN::FourierReconstructor::insert_slice(), EMAN::TestTomoImage::process_inplace(), EMAN::PhaseToMassCenterProcessor::process_inplace(), EMAN::ToMassCenterProcessor::process_inplace(), refalifn(), refalin3d_perturbquat(), EMAN::EMData::rotate_translate(), set_params(), set_params_inverse(), set_pre_trans(), EMAN::EMData::translate(), EMAN::RT3DSphereAligner::xform_align_nbest(), and EMAN::RT3DGridAligner::xform_align_nbest().

00876 {
00877         bool x_mirror = get_mirror();
00878 
00879         if (x_mirror) matrix[0][3] = -x;
00880         else matrix[0][3] = x;
00881         matrix[1][3] = y;
00882         matrix[2][3] = z;
00883 }

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.

Returns:
a Transforms, alt=54.73561, phi=45

Definition at line 83 of file transform.cpp.

References set_rotation(), and 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 177 of file transform.cpp.

References matrix.

Referenced by Transform().

00178 {
00179 //      transform_type = UNKNOWN;
00180         for(int i=0; i<3; ++i) {
00181                 for(int j=0; j<4; ++j) {
00182                         if(i==j) {
00183                                 matrix[i][j] = 1;
00184                         }
00185                         else {
00186                                 matrix[i][j] = 0;
00187                         }
00188                 }
00189         }
00190 }

template<typename Type>
Vec3f EMAN::Transform::transform const Vec3< Type > &  v  )  const [inline]
 

Transform a 3D vector using the internal transformation matrix.

Parameters:
v a three dimensional vector to be transformed
Returns:
the transformed vector

Definition at line 396 of file transform.h.

References v, and EMAN::Vec3f.

00396                                                                           {
00397 //                              assert_consistent_type(THREED); // Transform does the assertion
00398                                 return transform(v[0],v[1],v[2]);
00399                         }

Vec3f EMAN::Transform::transform const float &  x,
const float &  y,
const float &  z
const [inline]
 

Transform 3D coordinates using the internal transformation matrix.

Parameters:
x the x coordinate of the transformed point
y the y coordinate of the transformed point
z the z coordinate of the transformed point
Returns:
the transformed vector

Definition at line 382 of file transform.h.

References matrix(), EMAN::Vec3f, x, and y.

00382                                                                                                      {
00383 //                              assert_consistent_type(THREED);
00384                                 Vec3f ret;
00385                                 ret[0] = matrix[0][0] * x + matrix[0][1] * y + matrix[0][2] * z + matrix[0][3];
00386                                 ret[1] = matrix[1][0] * x + matrix[1][1] * y + matrix[1][2] * z + matrix[1][3];
00387                                 ret[2] = matrix[2][0] * x + matrix[2][1] * y + matrix[2][2] * z + matrix[2][3];
00388                                 return ret;
00389                         }

template<typename Type>
Vec2f EMAN::Transform::transform const Vec2< Type > &  v  )  const [inline]
 

Transform a 2D vector using the internal transformation matrix.

Parameters:
v a two dimensional vector to be transformed
Returns:
the transformed vector

Definition at line 372 of file transform.h.

References v, and EMAN::Vec2f.

00372                                                                           {
00373                                 return transform(v[0],v[1]);
00374                         }

Vec2f EMAN::Transform::transform const float &  x,
const float &  y
const [inline]
 

Transform 2D coordinates using the internal transformation matrix.

Parameters:
x the x coordinate of the transformed point
y the y coordinate of the transformed point
Returns:
the transformed vector

Definition at line 359 of file transform.h.

References matrix(), EMAN::Vec2f, x, and y.

Referenced by EMAN::EMData::get_rotated_clip(), and EMAN::operator *().

00359                                                                                      {
00360 //                              assert_valid_2d();
00361                                 Vec2f ret;
00362                                 ret[0] = matrix[0][0]*x + matrix[0][1]*y + matrix[0][3];
00363                                 ret[1] = matrix[1][0]*x + matrix[1][1]*y + matrix[1][3];
00364                                 return ret;
00365                         }

Transform Transform::transpose  )  const
 

Get the transpose of this transformation matrix.

Returns:
the transpose of this transformation matrix

Definition at line 1152 of file transform.cpp.

References t, and transpose_inplace().

Referenced by EMAN::PointArray::right_transform(), and EMAN::PDBReader::right_transform().

01152                                      {
01153         Transform t(*this);
01154         t.transpose_inplace();
01155         return t;
01156 }

void Transform::transpose_inplace  ) 
 

Get the transpose of this transformation matrix.

Definition at line 1139 of file transform.cpp.

References matrix.

Referenced by transpose().

01139                                   {
01140         float tempij;
01141         for (int i = 0; i < 3; i++) {
01142                 for (int j = 0; j < i; j++) {
01143                         if (i != j) {
01144                                 tempij= matrix[i][j];
01145                                 matrix[i][j] = matrix[j][i];
01146                                 matrix[j][i] = tempij;
01147                         }
01148                 }
01149         }
01150 }


Member Data Documentation

const float Transform::ERR_LIMIT = 0.000001f [static]
 

Definition at line 61 of file transform.cpp.

Referenced by get_determinant(), get_scale(), get_scale_and_mirror(), get_trans(), is_identity(), and set_scale().

float EMAN::Transform::matrix[3][4] [private]
 

Definition at line 442 of file transform.h.

Referenced by assert_valid_2d(), copy_matrix_into_array(), get_determinant(), get_matrix(), get_rotation(), get_trans(), get_trans_2d(), invert(), is_identity(), operator=(), operator==(), orthogonalize(), set_matrix(), set_mirror(), set_rotation(), set_scale(), set_trans(), to_identity(), 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 63 of file transform.cpp.

Referenced by detect_problem_keys(), and init_permissable_keys().

vector< string > Transform::permissable_3d_not_rot [static, private]
 

Definition at line 64 of file transform.cpp.

Referenced by detect_problem_keys(), and init_permissable_keys().

map< string, vector< string > > Transform::permissable_rot_keys [static, private]
 

Definition at line 65 of file transform.cpp.

Referenced by detect_problem_keys(), and init_permissable_keys().


The documentation for this class was generated from the following files:
Generated on Thu Mar 10 23:01:06 2011 for EMAN2 by  doxygen 1.3.9.1