AdmmKernel Class Reference

#include <reconstruct_ADMM.h>

Collaboration diagram for AdmmKernel:

Public Member Functions
void	initializeKernel (double alpha, double a, double _projectionStep)
double	projectionValueAt (double u, double v)
void	convolveKernelWithItself (double _autocorrStep)
void	applyCTFToKernelAutocorrelation (CTFDescription &ctf, double Ts, MultidimArray< double > &autocorrelationWithCTF)
void	getKernelAutocorrelation (MultidimArray< double > &autocorrelation)
void	computeGradient (MultidimArray< double > &gradient, char direction, bool adjoint=false)
void	computeKernel3D (MultidimArray< double > &kernel)

Public Attributes
double	supp
double	projectionStep
double	autocorrStep
double	alpha
Matrix1D< double >	projectionProfile
MultidimArray< std::complex< double > >	FourierProjectionAutocorr
FourierTransformer	transformer
MultidimArray< double >	projectionAutocorrWithCTF

Detailed Description

Definition at line 40 of file reconstruct_ADMM.h.

Member Function Documentation

applyCTFToKernelAutocorrelation()

void AdmmKernel::applyCTFToKernelAutocorrelation	(	CTFDescription &	ctf,
		double	Ts,
		MultidimArray< double > &	autocorrelationWithCTF
	)

Definition at line 764 of file reconstruct_ADMM.cpp.

{
    double dig2cont=1.0/Ts;
    double wx, wy;
    auto xdim=(int)XSIZE(projectionAutocorrWithCTF);
    int xdim_2=xdim/2;
    double ixdim=1.0/xdim;
    double maxFreq=2*autocorrStep;
    double maxFreq2=maxFreq*maxFreq;
    // Initialize Fourier transform to 0
    memset(&A2D_ELEM(transformer.fFourier,0,0),0,MULTIDIM_SIZE(transformer.fFourier)*sizeof(std::complex<double>));
    for (int i=((transformer.fFourier).yinit); i<=((transformer.fFourier).yinit + (int)(transformer.fFourier).ydim - 1); ++i)
    {
        FFT_IDX2DIGFREQ_FAST(i,xdim,xdim_2,ixdim,wy);
        if (fabs(wy)>maxFreq)
            continue;
        double wy2=wy*wy;
        wy*=dig2cont;
        for (int j=((transformer.fFourier).xinit); j<=((transformer.fFourier).xinit + (int)(transformer.fFourier).xdim - 1); ++j)
        {
            FFT_IDX2DIGFREQ_FAST(j,xdim,xdim_2,ixdim,wx);
            if (fabs(wx)>maxFreq)
                continue;
            double wx2=wx*wx;
            if (wy2+wx2>maxFreq2)
                continue;
            wx*=dig2cont;
            ctf.precomputeValues(wx,wy);
            // A2D_ELEM(transformer.fFourier,i,j)=A2D_ELEM(FourierProjectionAutocorr,i,j)*fabs(ctf.getValueAt());
            A2D_ELEM(transformer.fFourier,i,j)=fabs(ctf.getValueAt());
        }
    }
    transformer.inverseFourierTransform();
    autocorrelationWithCTF=projectionAutocorrWithCTF;
    CenterFFT(autocorrelationWithCTF, false);

// Debugging code
//    Image<double> save;
//    save()=autocorrelationWithCTF;
//    save.write("PPPautocorrelationWithCTF.xmp");
//    std::cout << "Press any key" << std::endl;
//    char c; std::cin>> c;
}

computeGradient()

void AdmmKernel::computeGradient	(	MultidimArray< double > &	gradient,
		char	direction,
		bool	adjoint = false
	)

Definition at line 808 of file reconstruct_ADMM.cpp.

{
    double supp2=supp*supp;
    double K=-1/(supp2*bessi2(alpha));
    gradient.initZeros();
    for (int k=-supp; k<=supp; ++k)
    {
        int k2=k*k;
        for (int i=-supp; i<=supp; ++i)
        {
            int i2=i*i;
            for (int j=-supp; j<=supp; ++j)
            {
                int j2=j*j;
                double r2=k2+i2+j2;
                if (r2>supp2)
                    continue;
                double z=alpha*sqrt(1.0-r2/supp2);
                double value=K*z*bessi1(z);
                switch (direction)
                {
                case 'x': value*=j; break;
                case 'y': value*=i; break;
                case 'z': value*=k; break;
                }
                if (adjoint)
                    value*=-1;
                A3D_ELEM(gradient,k,i,j)=value;
            }
        }
    }
}

computeKernel3D()

void AdmmKernel::computeKernel3D

(

MultidimArray< double > &

kernel

)

Compute the kernel in 3D space

Definition at line 841 of file reconstruct_ADMM.cpp.

{
    double supp2=supp*supp;
    double K=1/(bessi2(alpha)*alpha*alpha);
    kernel.initZeros();
    for (int k=-supp; k<=supp; ++k)
    {
        int k2=k*k;
        for (int i=-supp; i<=supp; ++i)
        {
            int i2=i*i;
            for (int j=-supp; j<=supp; ++j)
            {
                int j2=j*j;
                double r2=k2+i2+j2;
                if (r2>supp2)
                    continue;
                double z=alpha*sqrt(1.0-r2/supp2);
                double value=K*z*z*bessi2(z);
                A3D_ELEM(kernel,k,i,j)=value;
            }
        }
    }
}

convolveKernelWithItself()

void AdmmKernel::convolveKernelWithItself

(

double

_autocorrStep

)

Definition at line 723 of file reconstruct_ADMM.cpp.

{
    autocorrStep=_autocorrStep;
    size_t length=ceil(2.5*supp/autocorrStep)+1;
    projectionAutocorrWithCTF.initZeros(2*length+1,2*length+1);
    projectionAutocorrWithCTF.setXmippOrigin();
    double isupp=1.0/supp;
    double K=supp/bessi2(alpha)*sqrt(2.0*PI/alpha);
    FOR_ALL_ELEMENTS_IN_ARRAY2D(projectionAutocorrWithCTF)
    {
        double s=sqrt(i*i+j*j)*autocorrStep;
        double tmp=s*isupp;
        if (tmp<1)
        {
            tmp=sqrt(1.0 - tmp*tmp);
            A2D_ELEM(projectionAutocorrWithCTF,i,j)=K*pow(tmp,2.5)*bessi2_5(alpha*tmp);
        }
    }

    transformer.FourierTransform(projectionAutocorrWithCTF,FourierProjectionAutocorr);
    K=MULTIDIM_SIZE(projectionAutocorrWithCTF)*autocorrStep*autocorrStep;
    FOR_ALL_ELEMENTS_IN_ARRAY2D(FourierProjectionAutocorr)
        A2D_ELEM(FourierProjectionAutocorr,i,j)*=A2D_ELEM(FourierProjectionAutocorr,i,j)*K;

//  Debugging code
//    transformer.inverseFourierTransform();
//    Image<double> save;
//    save()=projectionAutocorrWithCTF;
//    save.write("PPPautocorrelationWithoutCTF.xmp");
//    std::cout << "Press any key" << std::endl;
//    char c; std::cin>> c;
}

getKernelAutocorrelation()

void AdmmKernel::getKernelAutocorrelation

(

MultidimArray< double > &

autocorrelation

)

Definition at line 756 of file reconstruct_ADMM.cpp.

{
    transformer.fFourier=FourierProjectionAutocorr;
    transformer.inverseFourierTransform();
    autocorrelation=projectionAutocorrWithCTF;
    CenterFFT(autocorrelation, false);
}

initializeKernel()

void AdmmKernel::initializeKernel	(	double	alpha,
		double	a,
		double	_projectionStep
	)

Definition at line 688 of file reconstruct_ADMM.cpp.

{
    projectionStep=_projectionStep;
    supp=a;
    alpha=_alpha;
    size_t length=ceil(a/projectionStep)+1;
    projectionProfile.initZeros(length);
    double ia=1.0/a;
    double K=a/bessi2(alpha)*sqrt(2.0*PI/alpha);
    FOR_ALL_ELEMENTS_IN_MATRIX1D(projectionProfile)
    {
        double s=i*projectionStep;
        double tmp=s*ia;
        tmp=sqrt(1.0 - tmp*tmp);
        VEC_ELEM(projectionProfile,i)=K*pow(tmp,2.5)*bessi2_5(alpha*tmp);
        //  function p = KaiserBesselProjection(m, alpha, a, s)
        //      tmp = sqrt(1 - (s/a).^2);
        //      p = a ./ besseli(m, alpha) .* sqrt(2*pi/alpha) .* tmp.^(m+0.5) .* besseli(m+0.5, alpha*tmp);
        //  end
    }
}

projectionValueAt()

double AdmmKernel::projectionValueAt	(	double	u,
		double	v
	)

Definition at line 710 of file reconstruct_ADMM.cpp.

                                                       {
    double r = sqrt(u*u+v*v);
    if (r > supp) {
        return 0.;
    } else {
        r = r/projectionStep ;
        int rmin = std::floor(r);
        int rmax = rmin+1    ;
        double p = r-rmin;
        return p * VEC_ELEM(projectionProfile,rmin) + (1-p) * VEC_ELEM(projectionProfile,rmax);
    }
}

Member Data Documentation

alpha

double AdmmKernel::alpha

Definition at line 46 of file reconstruct_ADMM.h.

autocorrStep

double AdmmKernel::autocorrStep

Definition at line 45 of file reconstruct_ADMM.h.

FourierProjectionAutocorr

MultidimArray< std::complex<double> > AdmmKernel::FourierProjectionAutocorr

Definition at line 48 of file reconstruct_ADMM.h.

projectionAutocorrWithCTF

MultidimArray<double> AdmmKernel::projectionAutocorrWithCTF

Definition at line 50 of file reconstruct_ADMM.h.

projectionProfile

Matrix1D<double> AdmmKernel::projectionProfile

Definition at line 47 of file reconstruct_ADMM.h.

projectionStep

double AdmmKernel::projectionStep

Definition at line 44 of file reconstruct_ADMM.h.

supp

double AdmmKernel::supp

Definition at line 43 of file reconstruct_ADMM.h.

transformer

FourierTransformer AdmmKernel::transformer

Definition at line 49 of file reconstruct_ADMM.h.

---

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

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