EMAN2
Public Member Functions | Static Public Member Functions | Static Public Attributes | Private Member Functions | Private Attributes | Static Private Attributes
EMAN::Transform Class Reference

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

#include <transform.h>

Collaboration diagram for EMAN::Transform:
Collaboration graph
[legend]

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.
void rotate_origin (const Transform &by)
 Increment the rotation by multipling the rotation bit of the argument transfrom by the rotation part of the current transfrom.
void rotate_origin_newBasis (const Transform &tcs, const float &omega, const float &n1, const float &n2, const float &n3)
 Increment the rotation by multipling the rotation bit of the argument transfrom by the rotation part of the current transfrom This version rotates in the standard coordinate system, even after it have been modified by tcs.
void rotate (const Transform &by)
 Increment the rotation by multipling the rotation bit of the argument transfrom by the current transfrom.
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.
void translate (const float &tx, const float &ty, const float &tz=0)
 Increment the current translation by tx, ty, tz.
void translate (const Vec3f &v)
 Increment the current translation using vec3f& v.
void translate (const Vec2f &v)
 Increment the current translation using vec2f& v.
void translate_newBasis (const Transform &tcs, const float &tx, const float &ty, const float &tz=0)
 Increment the current translation by tx, ty, tz using a non standard basis Actualy what it does is remove the effect of tcs when a composite transfrom tcs*t (where t is the current transform) This function is used in the scenegraph.
void translate (const Transform &tcs, const Vec3f &v)
 Increment the current translation using vec3f& v and a non standard basis.
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.
void scale (const float &scale)
 Increment the scale.
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.
bool is_rot_identity () const
 Returns whethers or this matrix rotation 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.
vector< float > get_matrix_4x4 () const
 Get the 4x4 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.
vector< Transformget_sym_proj (const string &sym) const
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

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


Constructor & Destructor Documentation

Transform::Transform ( )

Default constructor Internal matrix is the identity.

Definition at line 104 of file transform.cpp.

References to_identity().

Referenced by rotate_origin_newBasis(), set_params_inverse(), and translate_newBasis().

{
        to_identity();
}
Transform::Transform ( const Transform rhs)

Copy constructor.

Parameters:
rhsthe object to be copied

Definition at line 109 of file transform.cpp.

{
        *this = that;
}
Transform::Transform ( const Dict d)

Construction using a dictionary.

Parameters:
dthe dictionary containing the parameters

Definition at line 135 of file transform.cpp.

References set_params(), and to_identity().

                                   {
        to_identity();
        set_params(d);
}
Transform::Transform ( const float  array[12])

Construction using an array of floats.

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

Definition at line 141 of file transform.cpp.

References matrix.

                                          {
        memcpy(matrix,array,12*sizeof(float));
}
Transform::Transform ( const vector< float >  array)

Construction using a vector of size 12.

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

Definition at line 145 of file transform.cpp.

References set_matrix().

{
        set_matrix(array);
}
EMAN::Transform::~Transform ( ) [inline]

Definition at line 131 of file transform.h.

{ }

Member Function Documentation

void Transform::assert_valid_2d ( ) const [private]

Definition at line 1319 of file transform.cpp.

References ERR_LIMIT, matrix, and UnexpectedBehaviorException.

Referenced by get_rotation(), and set_rotation().

                                      {
        int rotation_error = 0;
        int translation_error = 0;
        if (fabs(matrix[2][0]) > ERR_LIMIT) rotation_error++;
        if (fabs(matrix[2][1]) > ERR_LIMIT) rotation_error++;
        if (fabs(matrix[2][3]) > ERR_LIMIT) translation_error++;
        if (fabs(matrix[0][2]) > ERR_LIMIT) rotation_error++;
        if (fabs(matrix[1][2]) > ERR_LIMIT) rotation_error++;
//      if (fabs(matrix[2][2]-1.0) >ERR_LIMIT) rotation_error++; 
        if ( translation_error && rotation_error ) {
                throw UnexpectedBehaviorException("Error, the internal matrix contains 3D rotations and 3D translations. This object can not be considered 2D");
        } else if ( translation_error ) {
                throw UnexpectedBehaviorException("Error, the internal matrix contains a non zero z component for a 3D translation. This object can not be considered 2D");
        }
        else if ( rotation_error ) {
                throw UnexpectedBehaviorException("Error, the internal matrix contains 3D rotations and this object can not be considered 2D");
        }

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

References matrix.

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

{ return matrix[r][c]; }
void Transform::copy_matrix_into_array ( float * const  array) const

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

                                                               {

        int idx = 0;
        for(int i=0; i<3; ++i) {
                for(int j=0; j<4; ++j) {
                        array[idx] = matrix[i][j];
                        idx ++;
                }
        }
}
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:
dthe dictionary that was the function argument of the set_params, set_rotation or the set_params_inv function
Exceptions:
InvalidParameterExceptionif the dictionary is invalid in anyway

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

                                                 {
        if (permissable_rot_keys.size() == 0 ) {
                init_permissable_keys();
        }

        vector<string> verification;
        vector<string> problem_keys;
        bool is_2d = false;
        if (d.has_key_ci("type") ) {
                string type = Util::str_to_lower((string)d["type"]);
                bool problem = false;
                if (permissable_rot_keys.find(type) == permissable_rot_keys.end() ) {
                        problem_keys.push_back(type);
                        problem = true;
                }
                if ( !problem ) {
                        vector<string> perm = permissable_rot_keys[type];
                        std::copy(perm.begin(),perm.end(),back_inserter(verification));

                        if ( type == "2d" ) {
                                is_2d = true;
                                std::copy(permissable_2d_not_rot.begin(),permissable_2d_not_rot.end(),back_inserter(verification));
                        }
                }
        }
        if ( !is_2d ) {
                std::copy(permissable_3d_not_rot.begin(),permissable_3d_not_rot.end(),back_inserter(verification));
        }

        for (Dict::const_iterator it = d.begin(); it != d.end();  ++it) {
                if ( std::find(verification.begin(),verification.end(), it->first) == verification.end() ) {
                        problem_keys.push_back(it->first);
                }
        }

        if (problem_keys.size() != 0 ) {
                string error;
                if (problem_keys.size() == 1) {
                        error = "Transform Error: The \"" +problem_keys[0]+ "\" key is unsupported";
                } else {
                        error = "Transform Error: The ";
                        for(vector<string>::const_iterator cit = problem_keys.begin(); cit != problem_keys.end(); ++cit ) {
                                if ( cit != problem_keys.begin() ) {
                                        if (cit == (problem_keys.end() -1) ) error += " and ";
                                        else error += ", ";
                                }
                                error += "\"";
                                error += *cit;
                                error += "\"";
                        }
                        error += " keys are unsupported";
                }
                throw InvalidParameterException(error);
        }
}
float Transform::get_determinant ( ) const

Get the determinant of the matrix.

Returns:
the determinant

Definition at line 1232 of file transform.cpp.

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

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

{
        float det;
        double det2;
        det2  = matrix[0][0]*((double)matrix[1][1]*matrix[2][2]-(double)matrix[2][1]*matrix[1][2]);
        det2 -= matrix[0][1]*((double)matrix[1][0]*matrix[2][2]-(double)matrix[2][0]*matrix[1][2]);
        det2 += matrix[0][2]*((double)matrix[1][0]*matrix[2][1]-(double)matrix[2][0]*matrix[1][1]);

        det = (float)det2;
        Util::apply_precision(det,ERR_LIMIT);

        return det;
}
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 778 of file transform.cpp.

References get_rotation(), get_trans(), set_rotation(), and set_trans().

                                               {

        Dict rot = get_rotation("eman");
        rot["alt"] = 180.0f + static_cast<float>(rot["alt"]);
        rot["phi"] = 180.0f - static_cast<float>(rot["phi"]);
       // This is the same as new_alt= 180-alt, new_phi=-phi, new_az=180+az

        Transform ret(*this); // Is the identity
        ret.set_rotation(rot);

        Vec3f trans = get_trans();
        trans[0] = -trans[0];
        ret.set_trans(trans);

//      ret.set_mirror(self.get_mirror());

        return ret;
}
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 172 of file transform.cpp.

References matrix.

Referenced by EMAN::Util::BPCQ(), EMAN::EMData::extract_box(), rotate(), and rotate_origin().

{
        vector<float> ret(12);
        for(int i=0; i<3; ++i) {
                for(int j=0; j<4; ++j) {
                        ret[i*4+j] = matrix[i][j];
                }
        }
        return ret;
}
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:
ithe row number (starting at 0)
Returns:
the ith row

Definition at line 470 of file transform.h.

References matrix.

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

                                                                  {
                                return Vec3f(matrix[i][0], matrix[i][1], matrix[i][2]);
                        }
vector< float > Transform::get_matrix_4x4 ( ) const

Get the 4x4 transformation matrix using a vector.

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

Definition at line 183 of file transform.cpp.

References matrix.

{
        vector<float> ret(16);
        for(int i=0; i<3; ++i) {
                for(int j=0; j<4; ++j) {
                        ret[i*4+j] = matrix[i][j];
                }
        }
        ret[12] = 0.0; 
        ret[13] = 0.0;
        ret[14] = 0.0;
        ret[15] = 1.0;
        
        return ret;
}
bool Transform::get_mirror ( ) const

Query whether x_mirroring is occuring.

Returns:
whether x_mirroring is occuring

Definition at line 1203 of file transform.cpp.

References get_determinant().

Referenced by get_params(), get_params_inverse(), get_trans(), get_trans_2d(), EMAN::GaussFFTProjector::project3d(), set_mirror(), set_trans(), and translate().

                                 {
        float determinant = get_determinant();

        bool x_mirror = false;
        if ( determinant < 0 ) x_mirror = true;

        return x_mirror;

}
int Transform::get_nsym ( const string &  sym) [static]
Dict Transform::get_params ( const string &  euler_type) const

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

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

Definition at line 471 of file transform.cpp.

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

Referenced by EMAN::RefineAlignerCG::align(), EMAN::RefineAligner::align(), EMAN::Util::BPCQ(), compose_transform2(), diff_between_3D_parameters_angles(), EMAN::Util::get_transform_params(), EMAN::Util::multiref_polar_ali_2d_local(), EMAN::Util::multiref_polar_ali_2d_local_psi(), EMAN::Util::multiref_polar_ali_2d_peaklist_local(), EMAN::Util::multiref_polar_ali_helical_90_local(), EMAN::Util::multiref_polar_ali_helical_local(), EMAN::Util::multiref_polar_ali_helicon_90_local(), EMAN::Util::multiref_polar_ali_helicon_local(), EMAN::TomoTiltEdgeMaskProcessor::process_inplace(), and EMAN::SpiderIO::write_single_header().

                                                         {
        Dict params = get_rotation(euler_type);

        Vec3f v = get_trans();
        params["tx"] = v[0]; params["ty"] = v[1];

        string type = Util::str_to_lower(euler_type);
        if ( type != "2d") params["tz"] = v[2];

        float scale = get_scale();
        params["scale"] = scale;

        bool mirror = get_mirror();
        params["mirror"] = mirror;

        return params;
}
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_typethe euler type of the retrieved rotation
Returns:
a dictionary containing the parameters

Definition at line 491 of file transform.cpp.

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

                                                                 {
        Transform inv(inverse());

        Dict params = inv.get_rotation(euler_type);
        Vec3f v = inv.get_pre_trans();
        params["tx"] = v[0]; params["ty"] = v[1];

        string type = Util::str_to_lower(euler_type);
        if ( type != "2d") params["tz"] = v[2];

        float scale = inv.get_scale();
        params["scale"] = scale;

        bool mirror = inv.get_mirror();
        params["mirror"] = mirror;

        return params;
}
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 1053 of file transform.cpp.

References get_trans(), invert(), and set_trans().

Referenced by get_params_inverse().

{
        Transform T(*this);
        T.set_trans(0,0,0);
        T.invert();

        Transform soln  = T*(*this);
//      soln.printme();
        return soln.get_trans();
}
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 1064 of file transform.cpp.

References get_trans_2d(), invert(), and set_trans().

{
        Transform T(*this);
        T.set_trans(0,0,0);
        T.invert();

        Transform soln  = T*(*this);
//      soln.printme();
        return soln.get_trans_2d();
}
Dict Transform::get_rotation ( const string &  euler_type = "eman") const

Get a rotation in any Euler format.

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

Definition at line 815 of file transform.cpp.

References assert_valid_2d(), EMAN::EMConsts::deg2rad, get_scale_and_mirror(), InvalidStringException, matrix, phi, EMAN::EMConsts::rad2deg, scale(), 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_sym_proj(), get_vflip_transform(), EMAN::FourierReconstructorSimple2D::insert_slice(), main(), EMAN::ChaoProjector::project3d(), EMAN::FourierGriddingProjector::project3d(), set_pre_trans(), and EMAN::RT3DSymmetryAligner::xform_align_nbest().

{
        Dict result;

        //float max = 1 - ERR_LIMIT;
        float scale;
        bool x_mirror;
        get_scale_and_mirror(scale,x_mirror);
        if (scale == 0) throw UnexpectedBehaviorException("The determinant of the Transform is 0. This is unexpected.");

        double cosalt = matrix[2][2]/scale;
        double x_mirror_scale = (x_mirror ? -1.0f : 1.0f);
        double inv_scale = 1.0f/scale;

        double az  = 0;
        double alt = 0;
        double phi = 0;
        double phiS = 0;  // like az  (but in SPIDER ZYZ)
        double psiS = 0;  // like phi (but in SPIDER ZYZ)

        // get alt, az, phi in EMAN convention

        if (cosalt >= 1) {  // that is, alt close to 0
                        alt = 0;
                        az = 0;
                        phi = (double)EMConsts::rad2deg * atan2(x_mirror_scale*matrix[0][1], x_mirror_scale*matrix[0][0]);
        } else if (cosalt <= -1) {  // that is, alt close to 180
                        alt = 180;
                        az = 0;
                        phi = (double)EMConsts::rad2deg * atan2(-x_mirror_scale*matrix[0][1], x_mirror_scale*matrix[0][0]);
        } else {   // for non exceptional cases:  0 < alt < 180

                az  = (double)EMConsts::rad2deg * atan2(scale*matrix[2][0], -scale*matrix[2][1]);

                if (matrix[2][2]==0.0)
                        alt = 90.0;
                else
                        alt = (double)EMConsts::rad2deg * atan(sqrt((double)matrix[2][0]*matrix[2][0]+(double)matrix[2][1]*matrix[2][1])/fabs(matrix[2][2]));

                if (matrix[2][2] * scale < 0)
                        alt = 180.0f-alt;
                
                phi = (double)EMConsts::rad2deg * atan2(x_mirror_scale*(double)matrix[0][2], (double)matrix[1][2]);

        } // ends separate cases: alt close to 0, 180, or neither

        phi = phi-360.0*floor(phi/360.0);
        az  = az -360.0*floor(az/360.0);

//  get phiS, psiS (SPIDER)
        if (cosalt >= 1) {  // that is, alt close to 0
                phiS = 0;
                psiS = phi;
        } else if (cosalt <= -1) {  // that is, alt close to 180
                phiS = 0;
                psiS = phi + 180.0;
        } else {
                phiS = az  - 90.0;
                psiS = phi + 90.0;
        }

        phiS = phiS-360.0*floor(phiS/360.0);
        psiS = psiS-360.0*floor(psiS/360.0);

//   do some quaternionic stuff here
        double xtilt = 0;
        double ytilt = 0;
        double ztilt = 0;


        string type = Util::str_to_lower(euler_type);

        result["type"] = type;
        if (type == "2d") {
                assert_valid_2d();
                result["alpha"]  = phi;
        } else if (type == "eman") {
//              assert_consistent_type(THREED);
                result["az"]  = az;
                result["alt"] = alt;
                result["phi"] = phi;
        } else if (type == "imagic") {
//              assert_consistent_type(THREED);
                result["alpha"] = az;
                result["beta"]  = alt;
                result["gamma"] = phi;
        } else if (type == "spider") {
//              assert_consistent_type(THREED);
                result["phi"]   = phiS;  // The first Euler like az
                result["theta"] = alt;
                result["psi"]   = psiS;
        } else if (type == "mrc") {
//              assert_consistent_type(THREED);
                result["phi"]   = phiS;
                result["theta"] = alt;
                result["omega"] = psiS;
        } else if (type == "xyz") {               // need to double-check these 3 equations ********
//              assert_consistent_type(THREED);
                xtilt = atan2(-sin(EMConsts::deg2rad*phiS)*sin(EMConsts::deg2rad*alt),cos(EMConsts::deg2rad*alt));
                ytilt = asin(  cos(EMConsts::deg2rad*phiS)*sin(EMConsts::deg2rad*alt));
                ztilt = psiS*EMConsts::deg2rad - atan2(sin(xtilt), cos(xtilt) *sin(ytilt));

                xtilt *= EMConsts::rad2deg; ytilt *= EMConsts::rad2deg; ztilt *= EMConsts::rad2deg;
                xtilt = xtilt-360*.0*floor((xtilt+180.0)/360.0);
                ytilt = ytilt-360*.0*floor((ytilt+180.0)/360.0);  //already in range [-90,90] but anyway...
                ztilt = ztilt-360*.0*floor((ztilt+180.0)/360.0);

                result["xtilt"]  = xtilt;
                result["ytilt"]  = ytilt;
                result["ztilt"]  = ztilt;
        } else if ((type == "quaternion") || (type == "spin") ||  (type == "sgirot")) {
          
              // The cosOover2 is also e0
//              double nphi = (az-phi)/2.0;
//              double cosOover2 = cos((az+phi)*EMConsts::deg2rad/2.0) * cos(alt*EMConsts::deg2rad/2.0);
//              printf("%f %f %f",matrix[0][0],matrix[1][1],matrix[2][2]);
                double traceR = matrix[0][0]+matrix[1][1]+matrix[2][2]; // This should be 1 + 2 cos omega
                double cosomega =  (traceR-1.0)/2.0;
                if (cosomega>1.0) cosomega=1.0;
                if (cosomega<-1.0) cosomega=-1.0;
                
                  // matrix(x,y)-matrix(y,x) = 2 n_z   sin(omega) etc
                 // trace matrix = 1 + 2 cos(omega)
                double sinOover2= sqrt((1.0 -cosomega)/2.0);
                double cosOover2= sqrt(1.0 -sinOover2*sinOover2);
                double sinomega = 2* sinOover2*cosOover2; 
                double n1 = 0; double n2 = 0;   double n3 = 0;
                if (sinomega>0) {
                      n1 = (matrix[1][2]-matrix[2][1])/2.0/sinomega ;
                      n2 = (matrix[2][0]-matrix[0][2])/2.0/sinomega ;
                      n3 = (matrix[0][1]-matrix[1][0])/2.0/sinomega ;
                }
//              printf("traceR=%lf,OneMinusCosomega=%lf,sinOover2=%lf,cosOover2=%lf,sinomega=%lf,cosomega=%lf,n3=%lf \n",traceR,1-cosomega,sinOover2,cosOover2,sinomega,cosomega,n3);

                
                if (type == "quaternion"){
                    result["e0"] = cosOover2 ;
                    result["e1"] = sinOover2 * n1 ;
                    result["e2"] = sinOover2 * n2;
                    result["e3"] = sinOover2 * n3;
                }

                if (type == "spin"){
                    result["omega"] = EMConsts::rad2deg * acos(cosomega);
                    result["n1"] = n1;
                    result["n2"] = n2;
                    result["n3"] = n3;
                }

                if (type == "sgirot"){
                    result["q"] = EMConsts::rad2deg * acos(cosomega);
                    result["n1"] = n1;
                    result["n2"] = n2;
                    result["n3"] = n3;
                }
                    
        } else if (type == "matrix") {
//              assert_consistent_type(THREED);
                result["m11"] = x_mirror_scale*matrix[0][0]*inv_scale;
                result["m12"] = x_mirror_scale*matrix[0][1]*inv_scale;
                result["m13"] = x_mirror_scale*matrix[0][2]*inv_scale;
                result["m21"] = matrix[1][0]*inv_scale;
                result["m22"] = matrix[1][1]*inv_scale;
                result["m23"] = matrix[1][2]*inv_scale;
                result["m31"] = matrix[2][0]*inv_scale;
                result["m32"] = matrix[2][1]*inv_scale;
                result["m33"] = matrix[2][2]*inv_scale;
        } else {
                throw InvalidStringException(euler_type, "unknown Euler Type");
        }

        return result;
}
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 674 of file transform.cpp.

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

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

{
        Transform ret(*this);
        ret.set_scale(1.0);
        ret.set_mirror(false);
        ret.set_trans(0,0,0);
        //ret.orthogonalize(); // ?
        return ret;
}
float Transform::get_scale ( ) const

Get the scale that was applied.

Returns:
the scale factor

Definition at line 1098 of file transform.cpp.

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

Referenced by get_params(), get_params_inverse(), EMAN::TransformProcessor::process(), EMAN::TransformProcessor::process_inplace(), EMAN::GaussFFTProjector::project3d(), set_pre_trans(), and set_scale().

                                 {
        float determinant = get_determinant();
        if (determinant < 0 ) determinant *= -1;

        float scale = std::pow(determinant,1.0f/3.0f);
        int int_scale = static_cast<int>(scale);
        float scale_residual = scale-static_cast<float>(int_scale);
        if  ( scale_residual < ERR_LIMIT ) { scale = static_cast<float>(int_scale); };

        Util::apply_precision(scale, ERR_LIMIT);

        return scale;
}
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:
scalea reference to the value that will be assigned the scale value
x_mirrora reference to the value that will be assigned the x_mirror value

Definition at line 1213 of file transform.cpp.

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

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

                                                                       {

        float determinant = get_determinant();
        x_mirror = false;
        if ( determinant < 0 ) {
                x_mirror = true;
                determinant *= -1;
        }
        if (determinant != 1 ) {
                scale = std::pow(determinant,1.0f/3.0f);
                int int_scale = static_cast<int>(scale);
                float scale_residual = scale-static_cast<float>(int_scale);
                if  ( scale_residual < ERR_LIMIT ) { scale = static_cast<float>(int_scale); };
        }
        else scale = 1;

        Util::apply_precision(scale,ERR_LIMIT);
}
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 1341 of file transform.cpp.

References EMAN::Symmetry3D::get_sym().

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

{
        Symmetry3D* sym = Factory<Symmetry3D>::get(sym_name);
        Transform ret;
        ret = (*this) * sym->get_sym(n);
        delete sym;
        return ret;
}
vector< Transform > Transform::get_sym_proj ( const string &  sym) const

Definition at line 1350 of file transform.cpp.

References EMAN::EMConsts::deg2rad, EMAN::Symmetry3D::get_nsym(), get_rotation(), EMAN::Symmetry3D::get_sym(), matrix, set_rotation(), and t.

Referenced by EMAN::nn4_ctf_rectReconstructor::insert_padfft_slice(), EMAN::nn4_ctfReconstructor::insert_padfft_slice(), EMAN::nn4_rectReconstructor::insert_padfft_slice(), EMAN::nn4Reconstructor::insert_padfft_slice(), EMAN::nn4_ctfwReconstructor::insert_padfft_slice_weighted(), and EMAN::Util::shc().

{
        vector<Transform> ret;
        Transform t;
        Symmetry3D* sym = Factory<Symmetry3D>::get(sym_name);
        int nsym = sym->get_nsym();
        int n = nsym;
        
        if ((sym_name[0] == 'c' || sym_name[0] == 'd' ) &&  fabs(matrix[2][2]) < 1.e-6){
                
                Dict d1,d2;                             
                d2["theta"] = (double)90.0;
                d2["psi"] = (double)0.0;
                d2["phi"] = (double)0.0;
                d2["type"] = "spider";
                d1 = this->get_rotation("spider");
                
                if (sym_name[0] == 'c') {
                        if( nsym%2 == 0)        n = nsym/2;
                        
                        for (int k=0;k<n;k++) {                         
                                d2["phi"] = (double)d1["phi"] + k*double(360.0)/ nsym;
                                d2["psi"] = d1["psi"];
                                t.set_rotation(d2);
                                ret.push_back( t );
                        }
                                
                }
                else {
                        nsym = nsym/2;
                        
                        if (nsym%2 == 0) {
                                n = nsym;
                                float cos_phi = cos( EMConsts::deg2rad*360.0/2/nsym );
                        
                                for (int k=0;k<n;k++){
                                        
                                        if(k%2==0)      {
                                        
                                                d2["phi"] = (double)d1["phi"] + k/2*double(360.0)/ nsym;
                                                d2["psi"] = d1["psi"];
                                                t.set_rotation(d2);
                                                ret.push_back( t );     
                                        }
                                        else    {
                                                        
                                                if( ( fabs(1.0-matrix[2][0])>1.0e-6 )&& fabs( matrix[2][0]-cos_phi)>1.0e-6  ){
                                                        //cout<<"jumped into"<<endl;
                                                        d2["phi"] = k/2*double(360.0)/ nsym +180 - (double)d1["phi"];
                                                        d2["psi"] = (double)d1["psi"] + 180;
                                                        t.set_rotation(d2);
                                                        ret.push_back( t );
                                                }
                                        }
                                
                                }
                        }
                        
                        
                        
                        else    {
                                n = nsym*2;
                                float cos_phi = cos( EMConsts::deg2rad*360.0/4/nsym );
                                for (int k=0;k<n;k++){
                                        
                                        if(k%4==0)      {
                                        
                                                d2["phi"] = (double)d1["phi"] + k/4*360.0/ nsym;
                                                d2["psi"] = (double)d1["psi"];
                                                t.set_rotation(d2);
                                                ret.push_back( t );     
                                        }
                                        else if( k%4 ==1)       {
                                                if( ( fabs(1.0-matrix[2][0])>1.0e-6 )&& fabs( matrix[2][0]-cos_phi)>1.0e-6  ){
                                                
                                                        d2["phi"] = k/4*360.0/nsym + 360.0/2/nsym+180 - (double)d1["phi"];
                                                        d2["psi"] = (double)d1["psi"] + 180;
                                                        t.set_rotation(d2);
                                                        ret.push_back( t );
                                                }
                                
                                        }
                                        
                                        else if( k%4 ==2)       {
                                        
                                                d2["phi"] =  k/4*360.0/ nsym+360.0/2/nsym+180 + (double)d1["phi"];
                                                d2["psi"] = (double)d1["psi"];
                                                t.set_rotation(d2);
                                                ret.push_back( t );
                                
                                        }
                                        
                                        else if( k%4 ==3)       {
                                                if( ( fabs(1.0-matrix[2][0])>1.0e-6 )&& fabs( matrix[2][0]-cos_phi)>1.0e-6  ) {
                                                        d2["phi"] = k/4*360.0/nsym+ 2.0*360.0/2/nsym - (double)d1["phi"];
                                                        d2["psi"] = (double)d1["psi"] + 180;
                                                        t.set_rotation(d2);
                                                        ret.push_back( t );
                                                }
                                        }
                                
                                }
                        }
                        
                }
                
        }
        else {
                for (int k=0;k<nsym;k++) {
                        t =  sym->get_sym(k);
                        ret.push_back( (*this) * t );
                }
                
        }
        delete sym;
        return ret;
}
Vec3f Transform::get_trans ( ) const
Vec2f Transform::get_trans_2d ( ) const

Get the degenerant 2D post trans as a vec2f.

Returns:
the 2D translation vector

Definition at line 1041 of file transform.cpp.

References get_mirror(), matrix, and v.

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

{
        bool x_mirror = get_mirror();
        Vec2f v;
        if (x_mirror) v[0] = -matrix[0][3];
        else v[0] = matrix[0][3];
        v[1] = matrix[1][3];
        return v;
}
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 797 of file transform.cpp.

References get_rotation(), get_trans(), set_rotation(), and set_trans().

                                               {

        Dict rot = get_rotation("eman");
        rot["alt"] = 180.0f + static_cast<float>(rot["alt"]);
        rot["phi"] = - static_cast<float>(rot["phi"]);
        
       // This is the same as new_alt= 180-alt, new_phi=180-phi, new_az=180+az

        Transform ret(*this);
        ret.set_rotation(rot);

        Vec3f trans = get_trans();
        trans[1] = -trans[1];
        ret.set_trans(trans);

        return ret;
}
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

Doctor Steve says Phil's answer put the 2-fold in the wrong place based on the standard Virus convention (empirically). It was also 2 to 5 not 5 to 2. It is possible to rotate to a 2-fold directly from a 5-fold, though there are 2 possible orientations for the 2-2-2 convention, and this finds only one of them :^(

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 set_rotation(), and t.

                                 {
        Transform t;
        Dict d;
        d["type"] = "eman";
        d["phi"] = 0;
        d["az"] = 90.0f;
        d["alt"] = 31.717474; // 5 fold to the nearest 2-fold
        t.set_rotation(d);
        return t;
}
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 276 of file transform.cpp.

References permissable_2d_not_rot, permissable_3d_not_rot, and permissable_rot_keys.

Referenced by detect_problem_keys().

{

        permissable_2d_not_rot.push_back("tx");
        permissable_2d_not_rot.push_back("ty");
        permissable_2d_not_rot.push_back("scale");
        permissable_2d_not_rot.push_back("mirror");
        permissable_2d_not_rot.push_back("type");

        permissable_3d_not_rot.push_back("tx");
        permissable_3d_not_rot.push_back("ty");
        permissable_3d_not_rot.push_back("tz");
        permissable_3d_not_rot.push_back("scale");
        permissable_3d_not_rot.push_back("mirror");
        permissable_3d_not_rot.push_back("type");

        vector<string> tmp;
        tmp.push_back("alpha");
        permissable_rot_keys["2d"] = tmp;

        tmp.clear();
        tmp.push_back("alt");
        tmp.push_back("az");
        tmp.push_back("phi");
        permissable_rot_keys["eman"] = tmp;

        tmp.clear();
        tmp.push_back("psi");
        tmp.push_back("theta");
        tmp.push_back("phi");
        permissable_rot_keys["spider"] = tmp;

        tmp.clear();
        tmp.push_back("alpha");
        tmp.push_back("beta");
        tmp.push_back("gamma");
        permissable_rot_keys["imagic"] = tmp;

        tmp.clear();
        tmp.push_back("ztilt");
        tmp.push_back("xtilt");
        tmp.push_back("ytilt");
        permissable_rot_keys["xyz"] = tmp;

        tmp.clear();
        tmp.push_back("phi");
        tmp.push_back("theta");
        tmp.push_back("omega");
        permissable_rot_keys["mrc"] = tmp;

        tmp.clear();
        tmp.push_back("e0");
        tmp.push_back("e1");
        tmp.push_back("e2");
        tmp.push_back("e3");
        permissable_rot_keys["quaternion"] = tmp;

        tmp.clear();
        tmp.push_back("n1");
        tmp.push_back("n2");
        tmp.push_back("n3");
        tmp.push_back("omega");
        permissable_rot_keys["spin"] = tmp;

        tmp.clear();
        tmp.push_back("n1");
        tmp.push_back("n2");
        tmp.push_back("n3");
        tmp.push_back("q");
        permissable_rot_keys["sgirot"] = tmp;

        tmp.clear();
        tmp.push_back("m11");
        tmp.push_back("m12");
        tmp.push_back("m13");
        tmp.push_back("m21");
        tmp.push_back("m22");
        tmp.push_back("m23");
        tmp.push_back("m31");
        tmp.push_back("m32");
        tmp.push_back("m33");
        permissable_rot_keys["matrix"] = tmp;
}
Transform Transform::inverse ( ) const
void Transform::invert ( )

Get the inverse of this transformation matrix.

Returns:
the inverse of this transformation matrix

Definition at line 1246 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(), rotate_origin_newBasis(), set_params_inverse(), set_pre_trans(), translate_newBasis(), and EMAN::RT3DSphereAligner::xform_align_nbest().

                       {

        double m00 = matrix[0][0]; double m01=matrix[0][1]; double m02=matrix[0][2];
        double m10 = matrix[1][0]; double m11=matrix[1][1]; double m12=matrix[1][2];
        double m20 = matrix[2][0]; double m21=matrix[2][1]; double m22=matrix[2][2];
        double v0  = matrix[0][3]; double v1 =matrix[1][3]; double v2 =matrix[2][3];

        double cof00 = m11*m22-m12*m21;
        double cof11 = m22*m00-m20*m02;
        double cof22 = m00*m11-m01*m10;
        double cof01 = m10*m22-m20*m12;
        double cof02 = m10*m21-m20*m11;
        double cof12 = m00*m21-m01*m20;
        double cof10 = m01*m22-m02*m21;
        double cof20 = m01*m12-m02*m11;
        double cof21 = m00*m12-m10*m02;

        double det = m00* cof00 + m02* cof02 -m01*cof01;

        matrix[0][0] =   (float)(cof00/det);
        matrix[0][1] = - (float)(cof10/det);
        matrix[0][2] =   (float)(cof20/det);
        matrix[1][0] = - (float)(cof01/det);
        matrix[1][1] =   (float)(cof11/det);
        matrix[1][2] = - (float)(cof21/det);
        matrix[2][0] =   (float)(cof02/det);
        matrix[2][1] = - (float)(cof12/det);
        matrix[2][2] =   (float)(cof22/det);

        matrix[0][3] =  (float)((- cof00*v0 + cof10*v1 - cof20*v2)/det);
        matrix[1][3] =  (float)((  cof01*v0 - cof11*v1 + cof21*v2)/det);
        matrix[2][3] =  (float)((- cof02*v0 + cof12*v1 - cof22*v2)/det);
}
bool Transform::is_identity ( ) const

Returns whethers or this matrix is the identity.

Definition at line 213 of file transform.cpp.

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

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

                                  {
        for(int i=0; i<3; ++i) {
                for(int j=0; j<4; ++j) {
                        float c = matrix[i][j];
                        Util::apply_precision(c,ERR_LIMIT);
                        if(i==j) {
                                if (c != 1.0) return false;
                        }
                        else {
                                if (c != 0.0) return false;
                        }
                }
        }
        return true;
}
bool Transform::is_rot_identity ( ) const

Returns whethers or this matrix rotation is the identity.

Definition at line 229 of file transform.cpp.

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

                                      {
        for(int i=0; i<3; ++i) {
                for(int j=0; j<3; ++j) {
                        float c = matrix[i][j];
                        Util::apply_precision(c,ERR_LIMIT);
                        if(i==j) {
                                if (c != 1.0) return false;
                        }
                        else {
                                if (c != 0.0) return false;
                        }
                }
        }
        return true;
}
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 767 of file transform.cpp.

References set(), and t.

{
        Transform t(*this);
        for(unsigned int i = 0; i < 3; ++i) {
                for(unsigned int j = 0; j < 4; ++j) {
                        t.set(i,j,t[i][j]*-1);
                }
        }
        return t;
}
bool Transform::operator!= ( const Transform rhs) const

Unequality comparision operator.

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

Definition at line 131 of file transform.cpp.

References operator==().

                                                    {
        return !(operator==(rhs));
}
Transform & Transform::operator= ( const Transform that)

Assignment operator.

Parameters:
thatthat which this will become

Definition at line 114 of file transform.cpp.

References matrix.

                                                      {
        if (this != &that ) {
                memcpy(matrix,that.matrix,12*sizeof(float));
//              transform_type = that.transform_type;
        }
        return *this;
}
bool Transform::operator== ( const Transform rhs) const

Equality comparision operator.

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

Definition at line 122 of file transform.cpp.

References matrix.

Referenced by operator!=().

                                                    {
        if (memcmp(this->matrix, rhs.matrix, 3*4*sizeof(float)) == 0) {
                return true;
        }
        else {
                return false;
        }
}
float* EMAN::Transform::operator[] ( int  i) [inline]

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

Definition at line 410 of file transform.h.

References matrix.

{ return matrix[i]; }
const float* EMAN::Transform::operator[] ( int  i) const [inline]

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

Definition at line 415 of file transform.h.

References matrix.

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

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

{
        float scale;
        bool x_mirror;
        get_scale_and_mirror(scale,x_mirror);
        if (scale == 0) throw UnexpectedBehaviorException("The determinant of the Transform is 0. This is unexpected.");
        double inv_scale = 1.0/static_cast<double>(scale);
        double mirror_scale = (x_mirror == true ? -1.0:1.0);

        gsl_matrix * R = gsl_matrix_calloc(3,3);
        for ( unsigned int i = 0; i < 3; ++i )
        {
                for ( unsigned int j = 0; j < 3; ++j )
                {
                        if (i == 0 && mirror_scale != 1.0 ) {
                                gsl_matrix_set( R, i, j, static_cast<double>(matrix[i][j])*mirror_scale*inv_scale );
                        }
                        else {
                                gsl_matrix_set( R, i, j, static_cast<double>(matrix[i][j])*inv_scale );
                        }
                }
        }

        gsl_matrix * V = gsl_matrix_calloc(3,3);
        gsl_vector * S = gsl_vector_calloc(3);
        gsl_vector * work = gsl_vector_calloc(3);
        gsl_linalg_SV_decomp (R, V, S, work); // Now R is U of the SVD R = USV^T

        gsl_matrix * Soln = gsl_matrix_calloc(3,3);
        gsl_blas_dgemm (CblasNoTrans, CblasTrans, 1.0, R, V, 0.0, Soln);

        for ( unsigned int i = 0; i < 3; ++i )
        {
                for ( unsigned int j = 0; j < 3; ++j )
                {
                        matrix[i][j] = static_cast<float>( gsl_matrix_get(Soln,i,j) );
                }
        }

        // Apply scale if it existed previously
        if (scale != 1.0f) {
                for(int i=0; i<3; ++i) {
                        for(int j=0; j<3; ++j) {
                                matrix[i][j] *= scale;
                        }
                }
        }

        // Apply post x mirroring if it was applied previouslys
        if ( x_mirror ) {
                for(int j=0; j<3; ++j) {
                        matrix[0][j] *= -1.0f;
                }
        }

        gsl_matrix_free(V); gsl_matrix_free(R); gsl_matrix_free(Soln);
        gsl_vector_free(S); gsl_vector_free(work);
}
void EMAN::Transform::printme ( ) const [inline]

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

Definition at line 356 of file transform.h.

References matrix.

                                             {
                                printf("%8.6f %8.6f %8.6f %8.6f\n",matrix[0][0],matrix[0][1],matrix[0][2],matrix[0][3]);
                                printf("%8.6f %8.6f %8.6f %8.6f\n",matrix[1][0],matrix[1][1],matrix[1][2],matrix[1][3]);
                                printf("%8.6f %8.6f %8.6f %8.6f\n",matrix[2][0],matrix[2][1],matrix[2][2],matrix[2][3]);
                                printf("%8.6f %8.6f %8.6f %8.6f\n",0.0,0.0,0.0,1.0);

                        }
void Transform::rotate ( const Transform by)

Increment the rotation by multipling the rotation bit of the argument transfrom by the current transfrom.

Parameters:
rotationmultiplican, a tranform, R'', by which to multiply the current one, R' == R''R'

Definition at line 749 of file transform.cpp.

References get_matrix(), and matrix.

{
        vector<float> multmatrix = by.get_matrix();
        // First Multiply and put the result in a temp matrix
        Transform result;
        for (int i=0; i<3; i++) {
                for (int j=0; j<4; j++) {
                        result[i][j] = multmatrix[i*4]*matrix[0][j] +  multmatrix[i*4+1]*matrix[1][j] + multmatrix[i*4+2]*matrix[2][j];
                }
        }
        //Then put the result from the tmep matrix in the original one
        for (int i=0; i<3; i++) {
                for (int j=0; j<4; j++) {
                        matrix[i][j] = result[i][j];
                }
        }
}
void Transform::rotate_origin ( const Transform by)

Increment the rotation by multipling the rotation bit of the argument transfrom by the rotation part of the current transfrom.

Parameters:
bymultiplican, a tranform, R'', by which to multiply the current one, R' == R''R'

Definition at line 706 of file transform.cpp.

References get_matrix(), and matrix.

Referenced by rotate_origin_newBasis().

{
        vector<float> multmatrix = by.get_matrix();
        // First Multiply and put the result in a temp matrix
        Transform result;
        for (int i=0; i<3; i++) {
                for (int j=0; j<3; j++) {
                        result[i][j] = multmatrix[i*4]*matrix[0][j] +  multmatrix[i*4+1]*matrix[1][j] + multmatrix[i*4+2]*matrix[2][j];
                }
        }
        //Then put the result from the tmep matrix in the original one
        for (int i=0; i<3; i++) {
                for (int j=0; j<3; j++) {
                        matrix[i][j] = result[i][j];
                }
        }
}
void Transform::rotate_origin_newBasis ( const Transform tcs,
const float &  omega,
const float &  n1,
const float &  n2,
const float &  n3 
)

Increment the rotation by multipling the rotation bit of the argument transfrom by the rotation part of the current transfrom This version rotates in the standard coordinate system, even after it have been modified by tcs.

The effect is to undo the distortion casued by dcs. Useful in the scenegraph

Parameters:
tcs,thestansfrom that moves us to a non standard coordinate system
bymultiplican, a tranform, R'', by which to multiply the current one, R' == R''R'

Definition at line 724 of file transform.cpp.

References get_trans(), invert(), rotate_origin(), set_scale(), set_trans(), and Transform().

{
        //Get the rotational inverse
        Transform tcsinv = Transform(tcs);
        tcsinv.set_trans(0.0, 0.0, 0.0);
        tcsinv.set_scale(1.0);
        tcsinv.invert();
        
        //Get the current rotation
        Transform temp = Transform();
        temp.set_trans(n1, n2, n3);
        Transform cc = tcsinv*temp;
        Vec3f cctrans = cc.get_trans();
        
        //set the right rotation
        Dict spinrot = Dict();
        spinrot["type"] = "spin";
        spinrot["omega"] = omega;
        spinrot["n1"] = cctrans[0];
        spinrot["n2"] = cctrans[1];
        spinrot["n3"] = cctrans[2];
        Transform rightrot = Transform(spinrot);
        rotate_origin(rightrot);
}
void Transform::scale ( const float &  scale)

Increment the scale.

Parameters:
scalethe amount increment by

Definition at line 1112 of file transform.cpp.

References get_determinant(), and set_scale().

Referenced by get_params(), get_params_inverse(), get_rotation(), get_scale(), get_scale_and_mirror(), orthogonalize(), set_params(), set_params_inverse(), set_pre_trans(), and set_rotation().

{
        float determinant = get_determinant();
        if (determinant < 0) determinant *= -1.0f;
        float newscale = std::pow(determinant,1.0f/3.0f) + scale;
        if(newscale > 0.0001) set_scale(newscale); // If scale ~ 0 things blowup, so we need a little fudge factor
}
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 405 of file transform.h.

References matrix.

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

{ 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:
vthe transformation matrix stored as a vector - 3 rows of 4.

Definition at line 150 of file transform.cpp.

References InvalidParameterException, and matrix.

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

{
        if (v.size() != 12 ) throw InvalidParameterException("The construction array must be of size 12");

        for(int i=0; i<3; ++i) {
                for(int j=0; j<4; ++j) {
                        matrix[i][j] = v[i*4+j];
                }
        }
}
void Transform::set_mirror ( const bool  x_mirror)
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:
dthe dictionary containing the parameters

Definition at line 245 of file transform.cpp.

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

Referenced by Transform().

                                        {
        detect_problem_keys(d);

        if (d.has_key_ci("type") ) set_rotation(d);

        if (d.has_key_ci("scale")) {
                float scale = static_cast<float>(d.get_ci("scale"));
                set_scale(scale);
        }

        float dx=0,dy=0,dz=0;

        if (d.has_key_ci("tx")) dx = static_cast<float>(d.get_ci("tx"));
        if (d.has_key_ci("ty")) dy = static_cast<float>(d.get_ci("ty"));
        if (d.has_key_ci("tz")) dz = static_cast<float>(d.get_ci("tz"));

        if ( dx != 0.0 || dy != 0.0 || dz != 0.0 ) {
                set_trans(dx,dy,dz);
        }

        if (d.has_key_ci("mirror")) {
                EMObject e = d.get_ci("mirror");
                if ( (e.get_type() != EMObject::BOOL ) && (e.get_type() != EMObject::INT ) && (e.get_type() != EMObject::UNSIGNEDINT ) )
                        throw InvalidParameterException("Error, mirror must be a bool or an int");

                bool mirror = static_cast<bool>(e);
                set_mirror(mirror);
        }
}
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:
dthe dictionary containing the inverse parameters

Definition at line 417 of file transform.cpp.

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

                                                {
        detect_problem_keys(d);

        if (d.has_key_ci("type") ) set_rotation(d);

        float dx=0,dy=0,dz=0;
        if (d.has_key_ci("tx")) dx = static_cast<float>(d.get_ci("tx"));
        if (d.has_key_ci("ty")) dy = static_cast<float>(d.get_ci("ty"));
        if (d.has_key_ci("tz")) dz = static_cast<float>(d.get_ci("tz"));

        if ( (dx != 0.0 || dy != 0.0 || dz != 0.0) && d.has_key_ci("type") ) {
                Transform pre_trans;
                pre_trans.set_trans(dx,dy,dz);

                Transform tmp;
                tmp.set_rotation(d);

                if (d.has_key_ci("scale")) {
                        float scale = static_cast<float>(d.get_ci("scale"));
                        tmp.set_scale(scale);
                }

                Transform solution_trans = tmp*pre_trans;

                if (d.has_key_ci("scale")) {
                        Transform tmp;
                        float scale = static_cast<float>(d.get_ci("scale"));
                        tmp.set_scale(scale);
                        solution_trans = solution_trans*tmp;
                }

                tmp = Transform();
                tmp.set_rotation(d);
                solution_trans = solution_trans*tmp;
                set_trans(solution_trans.get_trans());
        }

        if (d.has_key_ci("scale")) {
                float scale = static_cast<float>(d.get_ci("scale"));
                set_scale(scale);
        }

        if (d.has_key_ci("mirror")) {
                EMObject e = d.get_ci("mirror");
                if ( (e.get_type() != EMObject::BOOL ) && (e.get_type() != EMObject::INT ) && (e.get_type() != EMObject::UNSIGNEDINT ) )
                        throw InvalidParameterException("Error, mirror must be a bool or an int");

                bool mirror = static_cast<bool>(e);
                set_mirror(mirror);
        }
        invert();
}
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:
vthe vector (Vec3f or Vec2f) that is the pre trans

Definition at line 559 of file transform.h.

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

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

                                                   {

                Transform tmp;
                Dict rot = get_rotation("eman");
                tmp.set_rotation(rot);

                float scale = get_scale();
                if (scale != 1.0 ) tmp.set_scale(scale);

                Transform trans;
                trans.set_trans(v);

                trans = tmp*trans;

                Transform tmp2;
                tmp2.set_rotation(rot);
                tmp2.invert(); // invert
                if (scale != 1.0 ) tmp2.set_scale(1.0f/scale);


                trans = trans*tmp2;

                set_trans(trans.get_trans());
        }
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:
rotationa dictionary containing all key-entry pair required of the associated Euler type

Definition at line 511 of file transform.cpp.

References assert_valid_2d(), EMAN::EMConsts::deg2rad, detect_problem_keys(), EMAN::Dict::get_ci(), get_scale_and_mirror(), EMAN::Dict::has_key_ci(), InvalidParameterException, InvalidStringException, matrix, phi, scale(), EMAN::Util::str_to_lower(), and UnexpectedBehaviorException.

Referenced by EMAN::SymAlignProcessor::align(), EMAN::RotateTranslateAligner::align(), EMAN::RotatePrecenterAligner::align(), get_hflip_transform(), get_sym_proj(), get_vflip_transform(), icos_5_to_2(), EMAN::FourierReconstructor::preprocess_slice(), EMAN::TestTomoImage::process_inplace(), EMAN::MrcIO::read_mrc_header(), EMAN::EMData::rotate_translate(), set_params(), set_params_inverse(), set_pre_trans(), set_rotation(), and tet_3_to_2().

{
        detect_problem_keys(rotation);
        string euler_type;

        if (!rotation.has_key_ci("type") ){
                        throw InvalidParameterException("argument dictionary does not contain the type key");
        }

        euler_type = static_cast<string>(rotation.get_ci("type"));// Warning, will throw


        double e0=0; double e1=0; double e2=0; double e3=0;
        double omega=0;
        double az  = 0;
        double alt = 0;
        double phi = 0;
        double cxtilt = 0;
        double sxtilt = 0;
        double cytilt = 0;
        double sytilt = 0;
        double cztilt = 0;
        double sztilt = 0;
        bool is_quaternion = 0;
        bool is_matrix = 0;
        bool is_xyz = 0;

        bool x_mirror;
        float scale;
        // Get these before anything changes so we can apply them again after the rotation is set
        get_scale_and_mirror(scale,x_mirror);
        if (scale == 0) throw UnexpectedBehaviorException("The determinant of the Transform is 0. This is unexpected.");

        string type = Util::str_to_lower(euler_type);
        if (type == "2d") {
                assert_valid_2d();
                az  = 0;
                alt = 0;
                phi = (double)rotation["alpha"] ;
        } else if ( type == "eman" ) {
//              validate_and_set_type(THREED);
                az  = (double)rotation["az"] ;
                alt = (double)rotation["alt"]  ;
                phi = (double)rotation["phi"] ;
        } else if ( type == "imagic" ) {
//              validate_and_set_type(THREED);
                az  = (double)rotation["alpha"] ;
                alt = (double)rotation["beta"]  ;
                phi = (double)rotation["gamma"] ;
        } else if ( type == "spider" ) {
//              validate_and_set_type(THREED);
                az =  (double)rotation["phi"]    + 90.0;
                alt = (double)rotation["theta"] ;
                phi = (double)rotation["psi"]    - 90.0;
        } else if ( type == "xyz" ) {
//              validate_and_set_type(THREED);
                is_xyz = 1;
                cxtilt = cos(EMConsts::deg2rad*(double)rotation["xtilt"]);
                sxtilt = sin(EMConsts::deg2rad*(double)rotation["xtilt"]);
                cytilt = cos(EMConsts::deg2rad*(double)rotation["ytilt"]);
                sytilt = sin(EMConsts::deg2rad*(double)rotation["ytilt"]);
                cztilt = cos(EMConsts::deg2rad*(double)rotation["ztilt"]);
                sztilt = sin(EMConsts::deg2rad*(double)rotation["ztilt"]);
        } else if ( type == "mrc" ) {
//              validate_and_set_type(THREED);
                az  = (double)rotation["phi"]   + 90.0f ;
                alt = (double)rotation["theta"] ;
                phi = (double)rotation["omega"] - 90.0f ;
        } else if ( type == "quaternion" ) {
//              validate_and_set_type(THREED);
                is_quaternion = 1;
                e0 = (double)rotation["e0"];
                e1 = (double)rotation["e1"];
                e2 = (double)rotation["e2"];
                e3 = (double)rotation["e3"];
        } else if ( type == "spin" ) {
//              validate_and_set_type(THREED);
                is_quaternion = 1;
                omega = (double)rotation["omega"];
                e0 = cos(omega*EMConsts::deg2rad/2.0);
                e1 = sin(omega*EMConsts::deg2rad/2.0) * (double)rotation["n1"];
                e2 = sin(omega*EMConsts::deg2rad/2.0) * (double)rotation["n2"];
                e3 = sin(omega*EMConsts::deg2rad/2.0) * (double)rotation["n3"];
        } else if ( type == "sgirot" ) {
//              validate_and_set_type(THREED);
                is_quaternion = 1;
                omega = (double)rotation["q"] ;
                e0 = cos(omega*EMConsts::deg2rad/2.0);
                e1 = sin(omega*EMConsts::deg2rad/2.0) * (double)rotation["n1"];
                e2 = sin(omega*EMConsts::deg2rad/2.0) * (double)rotation["n2"];
                e3 = sin(omega*EMConsts::deg2rad/2.0) * (double)rotation["n3"];
        } else if ( type == "matrix" ) {
                is_matrix = 1;
                matrix[0][0] = (float)rotation["m11"];
                matrix[0][1] = (float)rotation["m12"];
                matrix[0][2] = (float)rotation["m13"];
                matrix[1][0] = (float)rotation["m21"];
                matrix[1][1] = (float)rotation["m22"];
                matrix[1][2] = (float)rotation["m23"];
                matrix[2][0] = (float)rotation["m31"];
                matrix[2][1] = (float)rotation["m32"];
                matrix[2][2] = (float)rotation["m33"];
        } else {
//              transform_type = UNKNOWN;
                throw InvalidStringException(euler_type, "unknown Euler Type");
        }

        double azp  =  az*EMConsts::deg2rad;
        double altp = alt*EMConsts::deg2rad;
        double phip = phi*EMConsts::deg2rad;

        if (!is_quaternion && !is_matrix && !is_xyz) {
                matrix[0][0] =  (float)(cos(phip)*cos(azp) - cos(altp)*sin(azp)*sin(phip));
                matrix[0][1] =  (float)(cos(phip)*sin(azp) + cos(altp)*cos(azp)*sin(phip));
                matrix[0][2] =  (float)(sin(altp)*sin(phip));
                matrix[1][0] =  (float)(-sin(phip)*cos(azp) - cos(altp)*sin(azp)*cos(phip));
                matrix[1][1] =  (float)(-sin(phip)*sin(azp) + cos(altp)*cos(azp)*cos(phip));
                matrix[1][2] =  (float)(sin(altp)*cos(phip));
                matrix[2][0] =  (float)(sin(altp)*sin(azp));
                matrix[2][1] =  (float)(-sin(altp)*cos(azp));
                matrix[2][2] =  (float)cos(altp);
        }
        if (is_quaternion){
                matrix[0][0] = (float)(e0 * e0 + e1 * e1 - e2 * e2 - e3 * e3);
                matrix[0][1] = (float)(2.0f * (e1 * e2 + e0 * e3));
                matrix[0][2] = (float)(2.0f * (e1 * e3 - e0 * e2));
                matrix[1][0] = (float)(2.0f * (e2 * e1 - e0 * e3));
                matrix[1][1] = (float)(e0 * e0 - e1 * e1 + e2 * e2 - e3 * e3);
                matrix[1][2] = (float)(2.0f * (e2 * e3 + e0 * e1));
                matrix[2][0] = (float)(2.0f * (e3 * e1 + e0 * e2));
                matrix[2][1] = (float)(2.0f * (e3 * e2 - e0 * e1));
                matrix[2][2] = (float)(e0 * e0 - e1 * e1 - e2 * e2 + e3 * e3);
                // keep in mind matrix[0][2] is M13 gives an e0 e2 piece, etc
        }
        if (is_xyz){
                matrix[0][0] =  (float)(cytilt*cztilt);
                matrix[0][1] =  (float)(cxtilt*sztilt+sxtilt*sytilt*cztilt);
                matrix[0][2] =  (float)(sxtilt*sztilt-cxtilt*sytilt*cztilt);
                matrix[1][0] =  (float)(-cytilt*sztilt);
                matrix[1][1] =  (float)(cxtilt*cztilt-sxtilt*sytilt*sztilt);
                matrix[1][2] =  (float)(sxtilt*cztilt+cxtilt*sytilt*sztilt);
                matrix[2][0] =  (float)(sytilt);
                matrix[2][1] =  (float)(-sxtilt*cytilt);
                matrix[2][2] =  (float)(cxtilt*cytilt);
        }
        
        // Apply scale if it existed previously
        if (scale != 1.0f) {
                for(int i=0; i<3; ++i) {
                        for(int j=0; j<3; ++j) {
                                matrix[i][j] *= scale;
                        }
                }
        }

        // Apply post x mirroring if it was applied previously
        if ( x_mirror ) {
                for(int j=0; j<3; ++j) {
                        matrix[0][j] *= -1.0f;
                }
        }
}
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:
vthe direction you want to solve for

Definition at line 684 of file transform.cpp.

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

{
        if ( v[0] == 0 && v[1] == 0 && v[2] == 0 )
                throw UnexpectedBehaviorException("Can't set rotation for the null vector");

        Vec3f v1(v);
        v1.normalize();

        double theta = acos(v1[2]); // in radians
        double psi = atan2(v1[1],-v1[0]);

        Dict d;
        d["theta"] = (double)EMConsts::rad2deg*theta;
        d["psi"] = (double)EMConsts::rad2deg*psi;
        d["phi"] = (double)0.0;
        d["type"] = "spider";

        set_rotation(d);


}
void Transform::set_scale ( const float &  scale)

Set the scale.

Parameters:
scalethe amount to scale by

Definition at line 1076 of file transform.cpp.

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

Referenced by EMAN::RefineAlignerCG::align(), EMAN::RefineAligner::align(), EMAN::ScaleAligner::align(), EMAN::ScaleAlignerABS::align_using_base(), 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(), rotate_origin_newBasis(), scale(), EMAN::EMData::scale(), set_params(), set_params_inverse(), and set_pre_trans().

                                                {
        if (new_scale <= 0) {
                throw InvalidValueException(new_scale,"The scale factor in a Transform object must be positive and non zero");
        }
        // Transform = MTSR (Mirroring, Translation, Scaling, Rotate)
        // So changing the scale boils down to this....

        float old_scale = get_scale();

        float n_scale = new_scale;
        Util::apply_precision(n_scale,ERR_LIMIT);

        float corrected_scale = n_scale/old_scale;
        if ( corrected_scale != 1.0 ) {
                for(int i = 0; i < 3;  ++i ) {
                        for(int j = 0; j < 3; ++j ) {
                                matrix[i][j] *= corrected_scale;
                        }
                }
        }
}
void EMAN::Transform::set_trans ( const Vec3f v) [inline]

Set the post translation component using a Vec3f.

Parameters:
vthe 3D translation vector

Definition at line 224 of file transform.h.

References set_trans().

Referenced by set_trans().

{ set_trans(v[0],v[1],v[2]); }
void EMAN::Transform::set_trans ( const Vec2f v) [inline]

Set the post translation component using a Vec2f.

Parameters:
vthe 2D translation vector

Definition at line 229 of file transform.h.

References set_trans().

Referenced by set_trans().

{ set_trans(v[0],v[1]); }
void Transform::set_trans ( const float &  x,
const float &  y,
const float &  z = 0 
)
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

Doctor Phil says: AltAngle= acos(-1/3.0)*90/pi;

Definition at line 89 of file transform.cpp.

References set_rotation(), and t.

                                {
        Transform t;
        Dict d;
        d["type"] = "eman";
        d["phi"] = 45.0f;
        d["az"] = 0.0f;
        d["alt"] = 54.73561f; // 3 fold to a 2 fold
        t.set_rotation(d);
        return t;
}
void Transform::to_identity ( )

Force the internal matrix to become the identity.

Definition at line 198 of file transform.cpp.

References matrix.

Referenced by Transform().

{
//      transform_type = UNKNOWN;
        for(int i=0; i<3; ++i) {
                for(int j=0; j<4; ++j) {
                        if(i==j) {
                                matrix[i][j] = 1;
                        }
                        else {
                                matrix[i][j] = 0;
                        }
                }
        }
}
Vec3f EMAN::Transform::transform ( const float &  x,
const float &  y,
const float &  z 
) const [inline]

Transform 3D coordinates using the internal transformation matrix.

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

Definition at line 445 of file transform.h.

References matrix.

                                                                                                     {
//                              assert_consistent_type(THREED);
                                Vec3f ret;
                                ret[0] = matrix[0][0] * x + matrix[0][1] * y + matrix[0][2] * z + matrix[0][3];
                                ret[1] = matrix[1][0] * x + matrix[1][1] * y + matrix[1][2] * z + matrix[1][3];
                                ret[2] = matrix[2][0] * x + matrix[2][1] * y + matrix[2][2] * z + matrix[2][3];
                                return ret;
                        }
template<typename Type >
Vec3f EMAN::Transform::transform ( const Vec3< Type > &  v) const [inline]

Transform a 3D vector using the internal transformation matrix.

Parameters:
va three dimensional vector to be transformed
Returns:
the transformed vector

Definition at line 459 of file transform.h.

References transform().

                                                                          {
//                              assert_consistent_type(THREED); // Transform does the assertion
                                return transform(v[0],v[1],v[2]);
                        }
Vec2f EMAN::Transform::transform ( const float &  x,
const float &  y 
) const [inline]

Transform 2D coordinates using the internal transformation matrix.

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

Definition at line 422 of file transform.h.

References matrix.

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

                                                                                     {
//                              assert_valid_2d();
                                Vec2f ret;
                                ret[0] = matrix[0][0]*x + matrix[0][1]*y + matrix[0][3];
                                ret[1] = matrix[1][0]*x + matrix[1][1]*y + matrix[1][3];
                                return ret;
                        }
template<typename Type >
Vec2f EMAN::Transform::transform ( const Vec2< Type > &  v) const [inline]

Transform a 2D vector using the internal transformation matrix.

Parameters:
va two dimensional vector to be transformed
Returns:
the transformed vector

Definition at line 435 of file transform.h.

References transform().

                                                                          {
                                return transform(v[0],v[1]);
                        }
void Transform::translate ( const float &  tx,
const float &  ty,
const float &  tz = 0 
)

Increment the current translation by tx, ty, tz.

Parameters:
txthe x incrementation
tythe y incrementation
tzthe z incrementation

Definition at line 1016 of file transform.cpp.

References get_mirror(), and matrix.

Referenced by translate_newBasis().

{
        bool x_mirror = get_mirror();
        if (x_mirror) matrix[0][3] = -matrix[0][3] + tx;
        else matrix[0][3] = matrix[0][3] + tx;
        matrix[1][3] = matrix[1][3] + ty;
        matrix[2][3] = matrix[2][3] + tz;
}
void EMAN::Transform::translate ( const Vec2f v) [inline]

Increment the current translation using vec2f& v.

Parameters:
vthe translation vector

Definition at line 251 of file transform.h.

References translate().

Referenced by translate().

{ translate(v[0],v[1]); }
void EMAN::Transform::translate ( const Transform tcs,
const Vec3f v 
) [inline]

Increment the current translation using vec3f& v and a non standard basis.

Parameters:
vthe 3D translation vector

Definition at line 266 of file transform.h.

References translate_newBasis().

{ translate_newBasis(tcs, v[0],v[1],v[2]); }
void EMAN::Transform::translate ( const Vec3f v) [inline]

Increment the current translation using vec3f& v.

Parameters:
vthe 3D translation vector

Definition at line 246 of file transform.h.

References translate().

Referenced by translate().

{ translate(v[0],v[1],v[2]); }
void Transform::translate_newBasis ( const Transform tcs,
const float &  tx,
const float &  ty,
const float &  tz = 0 
)

Increment the current translation by tx, ty, tz using a non standard basis Actualy what it does is remove the effect of tcs when a composite transfrom tcs*t (where t is the current transform) This function is used in the scenegraph.

Parameters:
tcsthe transform specifing the new basis vectors
txthe x incrementation
tythe y incrementation
tzthe z incrementation

Definition at line 1025 of file transform.cpp.

References get_trans(), invert(), set_trans(), Transform(), and translate().

Referenced by translate().

{
        //Get the rotational inverse
        Transform tcsinv = Transform(tcs);
        tcsinv.set_trans(0.0, 0.0, 0.0);
        tcsinv.invert();
        
        //Now move the coordinate system
        Transform temp = Transform();
        temp.set_trans(tx, ty, tz);
        Transform nb_trans = tcsinv*temp;
        
        translate(nb_trans.get_trans());
        
}
Transform Transform::transpose ( ) const

Get the transpose of this transformation matrix.

Returns:
the transpose of this transformation matrix

Definition at line 1299 of file transform.cpp.

References t, and transpose_inplace().

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

                                     {
        Transform t(*this);
        t.transpose_inplace();
        return t;
}
void Transform::transpose_inplace ( )

Get the transpose of this transformation matrix.

Definition at line 1286 of file transform.cpp.

References matrix.

Referenced by transpose().

                                  {
        float tempij;
        for (int i = 0; i < 3; i++) {
                for (int j = 0; j < i; j++) {
                        if (i != j) {
                                tempij= matrix[i][j];
                                matrix[i][j] = matrix[j][i];
                                matrix[j][i] = tempij;
                        }
                }
        }
}

Member Data Documentation

const float Transform::ERR_LIMIT = 0.000001f [static]
float EMAN::Transform::matrix[3][4] [private]
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 511 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 512 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 513 of file transform.h.

Referenced by detect_problem_keys(), and init_permissable_keys().


The documentation for this class was generated from the following files: