Alignment Class Reference

#include <tomo_align_tilt_series.h>

Collaboration diagram for Alignment:

Classes
struct	AlignmentEstimation
class	AProgAlignSignificant
class	ARotationEstimator
class	AShiftCorrEstimator
class	AShiftEstimator
class	CudaRotPolarEstimator
class	CudaShiftCorrEstimator
class	IterativeAlignmentEstimator
class	PolarRotationEstimator
class	ProgAlignSignificantGPU
class	RotationEstimationSetting
class	ShiftCorrEstimator

Public Types
enum	AlignType { AlignType::None, AlignType::OneToN, AlignType::MToN, AlignType::Consecutive }

Public Member Functions
	Alignment (const ProgTomographAlignment *_prm)
void	clear ()
Alignment &	operator= (const Alignment &op)
double	optimizeGivenAxisDirection ()
void	computeGeometryDependentOfAxis ()
void	computeGeometryDependentOfRotation ()
double	computeError () const
void	computeErrorForLandmarks ()
void	updateModel ()

Public Attributes
const ProgTomographAlignment *	prm
std::vector< Matrix1D< double > >	di
std::vector< Matrix1D< double > >	rj
Matrix1D< double >	psi
double	rot
double	tilt
Matrix1D< double >	raxis
int	Nimg
int	Nlandmark
std::vector< Matrix2D< double > >	Ai
std::vector< Matrix2D< double > >	Ait
std::vector< Matrix2D< double > >	Aip
std::vector< Matrix2D< double > >	Aipt
std::vector< Matrix1D< double > >	barri
std::vector< Matrix1D< double > >	diaxis
std::vector< Matrix2D< double > >	B1i
std::vector< Matrix2D< double > >	B2i
Matrix2D< double >	Binvraxis
Matrix2D< double >	allLandmarksPredictedX
Matrix2D< double >	allLandmarksPredictedY
MultidimArray< double >	errorLandmark

Friends
std::ostream &	operator<< (std::ostream &out, Alignment &alignment)

Detailed Description

Definition at line 295 of file tomo_align_tilt_series.h.

Member Enumeration Documentation

AlignType

enum Alignment::AlignType

strong

Enumerator
None
OneToN
MToN
Consecutive

Definition at line 31 of file align_type.h.

{ None, OneToN, MToN, Consecutive };

Constructor & Destructor Documentation

Alignment()

Alignment::Alignment ( const ProgTomographAlignment *  _prm )

inline

Definition at line 307 of file tomo_align_tilt_series.h.

    {
        prm=_prm;
        Nimg=MAT_XSIZE(_prm->allLandmarksX);
        Nlandmark=MAT_YSIZE(_prm->allLandmarksX);
        clear();
    }

Member Function Documentation

clear()

void Alignment::clear ( )

inline

Clear

Definition at line 316 of file tomo_align_tilt_series.h.

    {
        psi.initZeros(Nimg);
        rot=90;
        tilt=90;
        raxis.initZeros(3);

        Ai.clear();
        Ait.clear();
        Aip.clear();
        Aipt.clear();
        di.clear();
        diaxis.clear();
        B1i.clear();
        B2i.clear();
        barri.clear();
        Matrix2D<double> dummy23(2,3);
        Matrix2D<double> dummy32(3,2);
        Matrix2D<double> dummy33(3,3);
        Matrix2D<double> dummy22(2,2);
        Matrix1D<double> dummy2(2);
        Matrix1D<double> dummy3(3);
        for (int i=0; i<Nimg; i++)
        {
            Ai.push_back(dummy23);
            Ait.push_back(dummy32);
            Aip.push_back(dummy23);
            Aipt.push_back(dummy32);
            di.push_back(dummy2);
            diaxis.push_back(dummy2);
            B1i.push_back(dummy33);
            B2i.push_back(dummy33);
            barri.push_back(dummy3);
        }
        allLandmarksPredictedX.initZeros(Nlandmark,Nimg);
        allLandmarksPredictedY.initZeros(Nlandmark,Nimg);
        errorLandmark.initZeros(Nlandmark);

        rj.clear();
        for (int j=0; j<Nlandmark; j++)
            rj.push_back(dummy3);
    }

computeError()

double Alignment::computeError

(

)

const

Compute error

Definition at line 2712 of file tomo_align_tilt_series.cpp.

{
    Matrix1D<double> pijp;
    double error=0;
    double N=0;
    for (int i=0; i<Nimg; i++)
    {
        int jjmax=prm->Vseti[i].size();
        for (int jj=0; jj<jjmax; jj++)
        {
            int j=prm->Vseti[i][jj];
            pijp=Ai[i]*rj[j]+di[i]+diaxis[i];
            MAT_ELEM(allLandmarksPredictedX,j,i)=XX(pijp);
            MAT_ELEM(allLandmarksPredictedY,j,i)=YY(pijp);
            double diffx=MAT_ELEM(prm->allLandmarksX,j,i)-XX(pijp);
            double diffy=MAT_ELEM(prm->allLandmarksY,j,i)-YY(pijp);
            error+=diffx*diffx+diffy*diffy;
            N++;
        }
    }
    return sqrt(error/N);
}

computeErrorForLandmarks()

void Alignment::computeErrorForLandmarks

(

)

Compute error for landmarks

Definition at line 2735 of file tomo_align_tilt_series.cpp.

{
    Nlandmark=MAT_YSIZE(prm->allLandmarksX);
    errorLandmark.initZeros(Nlandmark);
    for (int j=0; j<Nlandmark; j++)
    {
        int counterj=0;
        for (int i=0; i<Nimg; i++)
        {
            if (prm->allLandmarksX(j,i)!=XSIZE(*(prm->img[0])))
            {
                double diffx=MAT_ELEM(prm->allLandmarksX,j,i)-
                             MAT_ELEM(allLandmarksPredictedX,j,i);
                double diffy=MAT_ELEM(prm->allLandmarksY,j,i)-
                             MAT_ELEM(allLandmarksPredictedY,j,i);
                errorLandmark(j)+=sqrt(diffx*diffx+diffy*diffy);
                counterj++;
            }
        }
        errorLandmark(j)/=counterj;
    }
}

computeGeometryDependentOfAxis()

void Alignment::computeGeometryDependentOfAxis

(

)

Compute Aip, and its transpose

Definition at line 2647 of file tomo_align_tilt_series.cpp.

{
    Matrix1D<double> axis;
    Euler_direction(rot, tilt, 0, axis);

    // Compute Aip, Aipt
    Matrix2D<double> Raxis, I, Ip, B;
    Binvraxis.initZeros(3,3);
    I.initIdentity(3);
    Ip=I;
    Ip(2,2)=0;
    for (int i=0; i<Nimg; i++)
    {
        rotation3DMatrix(prm->tiltList[i],axis,Raxis,false);
        Aipt[i](0,0)=Aip[i](0,0)=Raxis(0,0);
        Aipt[i](1,0)=Aip[i](0,1)=Raxis(0,1);
        Aipt[i](2,0)=Aip[i](0,2)=Raxis(0,2);
        Aipt[i](0,1)=Aip[i](1,0)=Raxis(1,0);
        Aipt[i](1,1)=Aip[i](1,1)=Raxis(1,1);
        Aipt[i](2,1)=Aip[i](1,2)=Raxis(1,2);

        B=I-Raxis;
        B=B.transpose()*Ip*B;
        Binvraxis+=((double)prm->ni(i))*(B+B.transpose());

        B1i[i]=I-Raxis.transpose();
        B2i[i]=B1i[i]*Ip*Raxis;
        /*COSS:
                std::cout << "Dependent on axis i=" << i << std::endl
                          << "tilt[i]=" << prm->tiltList[i] << std::endl
                          << "Raxis=\n" << Raxis
                          << "Aip[i]=\n" << Aip[i]
                          << "B1i[i]=\n" << B1i[i] << std::endl
                          << "B2i[i]=\n" << B2i[i] << std::endl;
        */
    }
    Binvraxis=Binvraxis.inv();
}

computeGeometryDependentOfRotation()

void Alignment::computeGeometryDependentOfRotation

(

)

Compute Ai, and its transposes

Definition at line 2687 of file tomo_align_tilt_series.cpp.

{
    // Compute Ai, Ait
    Matrix2D<double> Rinplane, I(2,3);
    I(0,0)=I(1,1)=1;
    for (int i=0; i<Nimg; i++)
    {
        rotation3DMatrix(-psi(i),'Z',Rinplane,false);
        B1i[i]=B1i[i]*Rinplane.transpose();

        Rinplane.resize(2,2);
        Ai[i]=Rinplane*Aip[i];
        Ait[i]=Ai[i].transpose();
        diaxis[i]=Rinplane*(I-Aip[i])*raxis;
        /* COSS:
                std::cout << "Dependent on rotation i=" << i << std::endl
                          << "B1i[i]=\n" << B1i[i] << std::endl
                          << "I-Raxis=\n" << (I-Aip[i])
                          << "raxis=" << raxis.transpose() << std::endl;
                std::cout << "diaxis[" << i << "]=" << diaxis[i].transpose() << std::endl;
        */
    }
}

operator=()

Alignment& Alignment::operator= ( const Alignment &  op )

inline

Assignment operator

Definition at line 360 of file tomo_align_tilt_series.h.

    {
        if (this!=&op)
        {
            prm=op.prm;
            Nimg=op.Nimg;
            Nlandmark=op.Nlandmark;
            psi=op.psi;
            rot=op.rot;
            tilt=op.tilt;
            Ai=op.Ai;
            Ait=op.Ait;
            Aip=op.Aip;
            Aipt=op.Aipt;
            diaxis=op.diaxis;
            B1i=op.B1i;
            B2i=op.B2i;
            raxis=op.raxis;
            di=op.di;
            rj=op.rj;
            barri=op.barri;
            allLandmarksPredictedX=op.allLandmarksPredictedX;
            allLandmarksPredictedY=op.allLandmarksPredictedY;
            errorLandmark=op.errorLandmark;
        }
        return *this;
    }

optimizeGivenAxisDirection()

double Alignment::optimizeGivenAxisDirection

(

)

Optimize for rot

Definition at line 2610 of file tomo_align_tilt_series.cpp.

{
    double bestError=1e38;
    bool firstIteration=true, finish=false;
    int Niterations=0;
    computeGeometryDependentOfAxis();
    do
    {
        computeGeometryDependentOfRotation();
        double error=computeError();
#ifdef DEBUG

        std::cout << "it=" << Niterations << " Error= " << error
        << " raxis=" << raxis.transpose() << std::endl;
#endif

        updateModel();
        Niterations++;
        if (firstIteration)
        {
            bestError=error;
            firstIteration=false;
        }
        else
        {
            finish=((error>bestError) || (Niterations>1000) ||
                    (ABS(error-bestError)/bestError<0.001)) && Niterations>20;
            if (error<bestError)
                bestError=error;
        }
    }
    while (!finish);
    return bestError;
}

updateModel()

void Alignment::updateModel

(

)

Update 3D model

Definition at line 2761 of file tomo_align_tilt_series.cpp.

{
    Matrix1D<double> pij(2), piN(2);
    Matrix2D<double> A(3,3), AitAi;
    Matrix1D<double> b(3);

#ifdef DEBUG

    std::cout << "Step 0: error=" << computeError() << std::endl;
#endif

    // Update the 3D positions of the landmarks
    for (int j=0; j<Nlandmark; j++)
    {
        A.initZeros();
        b.initZeros();
        // Compute the part corresponding to Vj
        int iimax=prm->Vsetj[j].size();
        for (int ii=0; ii<iimax; ii++)
        {
            int i=prm->Vsetj[j][ii];
            XX(pij)=MAT_ELEM(prm->allLandmarksX,j,i);
            YY(pij)=MAT_ELEM(prm->allLandmarksY,j,i);
            A+=Ait[i]*Ai[i];
            b+=Ait[i]*(pij-(di[i]+diaxis[i]));
        }

        // Update rj[j]
        rj[j]=A.inv()*b;
    }
#ifdef DEBUG
    std::cout << "Step 1: error=" << computeError() << std::endl;
#endif

    // Compute the average landmarks seen in each image
    for (int i=0; i<Nimg; i++)
    {
        barri[i].initZeros();
        for (int jj=0; jj<prm->ni(i); jj++)
        {
            int j=prm->Vseti[i][jj];
            barri[i]+=rj[j];
        }
        barri[i]/=prm->ni(i);
    }

    // Update shifts
    if (prm->isCapillar)
    {
        // Update the raxis
        Matrix1D<double> tmp2(3), Pij(3);
        for (int i=0; i<Nimg; i++)
        {
            for (int jj=0; jj<prm->ni(i); jj++)
            {
                int j=prm->Vseti[i][jj];
                XX(Pij)=MAT_ELEM(prm->allLandmarksX,j,i);
                YY(Pij)=MAT_ELEM(prm->allLandmarksY,j,i);
                tmp2+=B1i[i]*Pij-B2i[i]*rj[j];
            }
        }
        tmp2*=2.0;
        tmp2=Binvraxis*tmp2;
        XX(raxis)=XX(tmp2);
        YY(raxis)=YY(tmp2);
        computeGeometryDependentOfRotation();
    }
    else
    {
        // Update the individual di
        for (int i=0; i<Nimg; i++)
            if (i!=prm->iMinTilt)
            {
                di[i] = prm->barpi[i]-Ai[i]*barri[i]-diaxis[i];

                if (di[i].isAnyNaN())
                    di[i].initZeros();
            }
    }
#ifdef DEBUG
    std::cout << "Step 2: error=" << computeError() << std::endl;
#endif

    // Update rotations
    if (prm->psiMax>0)
    {
        Matrix2D<double> Ri(2,2);
        Matrix2D<double> Aiprj(2,1), Aiprjt(1,2), dim(2,1), Pij(2,1);
        for (int i=0; i<Nimg; i++)
        {
            Ri.initZeros();
            dim.fromVector(di[i]+diaxis[i]);
            for (int jj=0; jj<prm->ni(i); jj++)
            {
                int j=prm->Vseti[i][jj];
                Aiprj.fromVector(Aip[i]*rj[j]);
                Aiprjt=Aiprj.transpose();

                MAT_ELEM(Pij,0,0)=MAT_ELEM(prm->allLandmarksX,j,i);
                MAT_ELEM(Pij,1,0)=MAT_ELEM(prm->allLandmarksY,j,i);
                Ri+=(dim-Pij)*Aiprjt;
            }
            psi(i)=CLIP(RAD2DEG(atan(((Ri(0,1)-Ri(1,0))/(Ri(0,0)+Ri(1,1))))),
                        -(prm->psiMax),prm->psiMax);
            Matrix2D<double> Rinplane;
            rotation3DMatrix(-psi(i),'Z',Rinplane,false);
        }
    }
#ifdef DEBUG
    std::cout << "Step 3: error=" << computeError() << std::endl;
#endif

    // Update the rotation dependent part
    computeGeometryDependentOfRotation();
}

Friends And Related Function Documentation

operator<<

std::ostream& operator<< ( std::ostream &  out, Alignment &  alignment  )

friend

Print an alignment

Definition at line 2879 of file tomo_align_tilt_series.cpp.

{
    out << "Alignment parameters ===========================================\n"
    << "rot=" << alignment.rot << " tilt=" << alignment.tilt << std::endl
    << "raxis=" << alignment.raxis.transpose() << std::endl;
    out << "Images ---------------------------------------------------------\n";
    for (int i=0; i<alignment.Nimg; i++)
        out << "Image " << i << " psi= " << alignment.psi(i)
        << " di= " << alignment.di[i].transpose()
        << " diaxis= " << alignment.diaxis[i].transpose()
        << " (" << alignment.prm->ni(i) << ")\n";
    out << "Landmarks ------------------------------------------------------\n";
    alignment.computeErrorForLandmarks();
    for (int j=0; j<alignment.Nlandmark; j++)
        out << "Landmark " << j << " rj= " << alignment.rj[j].transpose()
        << " " << alignment.errorLandmark(j) << std::endl;
    return out;
}

Member Data Documentation

Ai

std::vector< Matrix2D<double> > Alignment::Ai

Definition at line 418 of file tomo_align_tilt_series.h.

Aip

std::vector< Matrix2D<double> > Alignment::Aip

Definition at line 424 of file tomo_align_tilt_series.h.

Aipt

std::vector< Matrix2D<double> > Alignment::Aipt

Definition at line 427 of file tomo_align_tilt_series.h.

Ait

std::vector< Matrix2D<double> > Alignment::Ait

Definition at line 421 of file tomo_align_tilt_series.h.

allLandmarksPredictedX

Matrix2D<double> Alignment::allLandmarksPredictedX

Definition at line 445 of file tomo_align_tilt_series.h.

allLandmarksPredictedY

Matrix2D<double> Alignment::allLandmarksPredictedY

Definition at line 448 of file tomo_align_tilt_series.h.

B1i

std::vector< Matrix2D<double> > Alignment::B1i

Definition at line 436 of file tomo_align_tilt_series.h.

B2i

std::vector< Matrix2D<double> > Alignment::B2i

Definition at line 439 of file tomo_align_tilt_series.h.

barri

std::vector< Matrix1D<double> > Alignment::barri

Definition at line 430 of file tomo_align_tilt_series.h.

Binvraxis

Matrix2D<double> Alignment::Binvraxis

Definition at line 442 of file tomo_align_tilt_series.h.

di

std::vector< Matrix1D<double> > Alignment::di

Definition at line 299 of file tomo_align_tilt_series.h.

diaxis

std::vector< Matrix1D<double> > Alignment::diaxis

Definition at line 433 of file tomo_align_tilt_series.h.

errorLandmark

MultidimArray<double> Alignment::errorLandmark

Definition at line 451 of file tomo_align_tilt_series.h.

Nimg

int Alignment::Nimg

Definition at line 412 of file tomo_align_tilt_series.h.

Nlandmark

int Alignment::Nlandmark

Definition at line 415 of file tomo_align_tilt_series.h.

prm

const ProgTomographAlignment* Alignment::prm

Definition at line 298 of file tomo_align_tilt_series.h.

psi

Matrix1D<double> Alignment::psi

Definition at line 301 of file tomo_align_tilt_series.h.

raxis

Matrix1D<double> Alignment::raxis

Definition at line 304 of file tomo_align_tilt_series.h.

rj

std::vector< Matrix1D<double> > Alignment::rj

Definition at line 300 of file tomo_align_tilt_series.h.

rot

double Alignment::rot

Definition at line 302 of file tomo_align_tilt_series.h.

tilt

double Alignment::tilt

Definition at line 303 of file tomo_align_tilt_series.h.

---

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

• xmipp/legacy/libraries/reconstruction/tomo_align_tilt_series.h
• xmipp/legacy/libraries/reconstruction/tomo_align_tilt_series.cpp