Monogenic Class Reference

#include <monogenic_signal.h>

Public Member Functions
void	proteinRadiusVolumeAndShellStatistics (const MultidimArray< int > &mask, int &radius, long &vol)
void	findCliffValue (MultidimArray< double > &inputmap, int &radius, int &radiuslimit, MultidimArray< int > &mask, double rsmooth)
Matrix1D< double >	fourierFreqVector (size_t dimarrayFourier, size_t dimarrayReal)
MultidimArray< double >	fourierFreqs_3D (const MultidimArray< std::complex< double > > &myfftV, const MultidimArray< double > &inputVol, Matrix1D< double > &freq_fourier_x, Matrix1D< double > &freq_fourier_y, Matrix1D< double > &freq_fourier_z)
void	resolution2eval (int &count_res, double step, double &resolution, double &last_resolution, double &freq, double &freqL, int &last_fourier_idx, int &volsize, bool &continueIter, bool &breakIter, double &sampling, double &maxRes)
void	amplitudeMonoSig3D_LPF (const MultidimArray< std::complex< double > > &myfftV, FourierTransformer &transformer_inv, MultidimArray< std::complex< double > > &fftVRiesz, MultidimArray< std::complex< double > > &fftVRiesz_aux, MultidimArray< double > &VRiesz, double freq, double freqH, double freqL, MultidimArray< double > &iu, Matrix1D< double > &freq_fourier_x, Matrix1D< double > &freq_fourier_y, Matrix1D< double > &freq_fourier_z, MultidimArray< double > &amplitude, int count, FileName fnDebug)
void	statisticsInBinaryMask2 (const MultidimArray< double > &volS, const MultidimArray< double > &volN, MultidimArray< int > &mask, double &meanS, double &sdS2, double &meanN, double &sdN2, double &significance, double &thr95, double &NS, double &NN)
void	statisticsInOutBinaryMask2 (const MultidimArray< double > &volS, MultidimArray< int > &mask, double &meanS, double &sdS2, double &meanN, double &sdN2, double &significance, double &thr95, double &NS, double &NN)
void	setLocalResolutionHalfMaps (const MultidimArray< double > &amplitudeMS, MultidimArray< int > &pMask, MultidimArray< double > &plocalResolutionMap, double thresholdNoise, double resolution, double resolution_2)
void	setLocalResolutionMap (const MultidimArray< double > &amplitudeMS, MultidimArray< int > &pMask, MultidimArray< double > &plocalResolutionMap, double thresholdNoise, double resolution, double resolution_2)
void	monogenicAmplitude_3D_Fourier (const MultidimArray< std::complex< double > > &myfftV, MultidimArray< double > &iu, MultidimArray< double > &amplitude, int numberOfThreads)
void	addNoise (MultidimArray< double > &vol, double mean, double stddev)
MultidimArray< double >	createDataTest (size_t xdim, size_t ydim, size_t zdim, double wavelength)

Detailed Description

Definition at line 40 of file monogenic_signal.h.

Member Function Documentation

addNoise()

void Monogenic::addNoise	(	MultidimArray< double > &	vol,
		double	mean,
		double	stddev
	)

Definition at line 666 of file monogenic_signal.cpp.

{

    std::random_device rd;
    std::default_random_engine generator(rd());
    std::normal_distribution<double> dist(mean, stddev);

    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(vol)
        DIRECT_MULTIDIM_ELEM(vol,n) += dist(generator);
}

amplitudeMonoSig3D_LPF()

void Monogenic::amplitudeMonoSig3D_LPF	(	const MultidimArray< std::complex< double > > &	myfftV,
		FourierTransformer &	transformer_inv,
		MultidimArray< std::complex< double > > &	fftVRiesz,
		MultidimArray< std::complex< double > > &	fftVRiesz_aux,
		MultidimArray< double > &	VRiesz,
		double	freq,
		double	freqH,
		double	freqL,
		MultidimArray< double > &	iu,
		Matrix1D< double > &	freq_fourier_x,
		Matrix1D< double > &	freq_fourier_y,
		Matrix1D< double > &	freq_fourier_z,
		MultidimArray< double > &	amplitude,
		int	count,
		FileName	fnDebug
	)

Definition at line 283 of file monogenic_signal.cpp.

{
//FIXME: use atf.h
//  FourierTransformer transformer_inv;

    fftVRiesz.initZeros(myfftV);
    fftVRiesz_aux.initZeros(myfftV);
    std::complex<double> J(0,1);

    // Filter the input volume and add it to amplitude
    long n=0;
    double ideltal=PI/(freq-freqH);

    for(size_t k=0; k<ZSIZE(myfftV); ++k)
    {
        for(size_t i=0; i<YSIZE(myfftV); ++i)
        {
            for(size_t j=0; j<XSIZE(myfftV); ++j)
            {
                double iun=DIRECT_MULTIDIM_ELEM(iu,n);
                double un=1.0/iun;
                if (freqH<=un && un<=freq)
                {
                    DIRECT_MULTIDIM_ELEM(fftVRiesz, n) = DIRECT_MULTIDIM_ELEM(myfftV, n);
                    DIRECT_MULTIDIM_ELEM(fftVRiesz, n) *= 0.5*(1+cos((un-freq)*ideltal));//H;
                    DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n) = -J;
                    DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n) *= DIRECT_MULTIDIM_ELEM(fftVRiesz, n);
                    DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n) *= iun;
                } else if (un>freq)
                {
                    DIRECT_MULTIDIM_ELEM(fftVRiesz, n) = DIRECT_MULTIDIM_ELEM(myfftV, n);
                    DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n) = -J;
                    DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n) *= DIRECT_MULTIDIM_ELEM(fftVRiesz, n);
                    DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n) *= iun;
                }
                ++n;
            }
        }
    }

    transformer_inv.inverseFourierTransform(fftVRiesz, amplitude);
    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(amplitude)
        DIRECT_MULTIDIM_ELEM(amplitude,n) *= DIRECT_MULTIDIM_ELEM(amplitude,n);

    //TODO: create a macro with these kind of code
    // Calculate first component of Riesz vector
    double ux;
    n=0;
    for(size_t k=0; k<ZSIZE(myfftV); ++k)
    {
        for(size_t i=0; i<YSIZE(myfftV); ++i)
        {
            for(size_t j=0; j<XSIZE(myfftV); ++j)
            {
                ux = VEC_ELEM(freq_fourier_x,j);
                DIRECT_MULTIDIM_ELEM(fftVRiesz, n) = ux*DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n);
                ++n;
            }
        }
    }

    transformer_inv.inverseFourierTransform(fftVRiesz, VRiesz);
    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(amplitude)
        DIRECT_MULTIDIM_ELEM(amplitude,n)+=DIRECT_MULTIDIM_ELEM(VRiesz,n)*DIRECT_MULTIDIM_ELEM(VRiesz,n);


    // Calculate second and third component of Riesz vector
    n=0;
    double uy;
    double uz;
    for(size_t k=0; k<ZSIZE(myfftV); ++k)
    {
        uz = VEC_ELEM(freq_fourier_z,k);
        for(size_t i=0; i<YSIZE(myfftV); ++i)
        {
            uy = VEC_ELEM(freq_fourier_y,i);
            for(size_t j=0; j<XSIZE(myfftV); ++j)
            {
                DIRECT_MULTIDIM_ELEM(fftVRiesz, n) = uz*DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n);
                DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n) = uy*DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n);
                ++n;
            }
        }
    }

    transformer_inv.inverseFourierTransform(fftVRiesz, VRiesz);

    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(amplitude)
        DIRECT_MULTIDIM_ELEM(amplitude,n) += DIRECT_MULTIDIM_ELEM(VRiesz,n)*DIRECT_MULTIDIM_ELEM(VRiesz,n);

    transformer_inv.inverseFourierTransform(fftVRiesz_aux, VRiesz);

    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(amplitude)
    {
        DIRECT_MULTIDIM_ELEM(amplitude,n)+=DIRECT_MULTIDIM_ELEM(VRiesz,n)*DIRECT_MULTIDIM_ELEM(VRiesz,n);
        DIRECT_MULTIDIM_ELEM(amplitude,n)=sqrt(DIRECT_MULTIDIM_ELEM(amplitude,n));
    }

    transformer_inv.FourierTransform(amplitude, fftVRiesz, false);

    double raised_w = PI/(freqL-freq);
    n=0;
    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(fftVRiesz)
    {
        double un=1.0/DIRECT_MULTIDIM_ELEM(iu,n);
        if (freqL>=un && un>=freq)
        {
            DIRECT_MULTIDIM_ELEM(fftVRiesz,n) *= 0.5*(1 + cos(raised_w*(un-freq)));
        }
        else
        {
            if (un>freqL)
            {
                DIRECT_MULTIDIM_ELEM(fftVRiesz,n) = 0;
            }
        }
    }

    transformer_inv.inverseFourierTransform();
//
//  saveImg = amplitude;
//  FileName iternumber;
//  iternumber = formatString("_Amplitude_new_%i.vol", count);
//  saveImg.write(fnDebug+iternumber);
//  saveImg.clear();
}

createDataTest()

MultidimArray< double > Monogenic::createDataTest	(	size_t	xdim,
		size_t	ydim,
		size_t	zdim,
		double	wavelength
	)

Definition at line 681 of file monogenic_signal.cpp.

{
    int siz_z;
    int siz_y;
    int siz_x;
    double x;
    double y;
    double z;

    siz_z = xdim/2;
    siz_y = ydim/2;
    siz_x = xdim/2;
    MultidimArray<double> testmap;
    testmap.initZeros(zdim, ydim, xdim);

    long n=0;
    for(int k=0; k<zdim; ++k)
    {
        z = (k - siz_z);
        z= z*z;
        for(int i=0; i<ydim; ++i)
        {
            y = (i - siz_y);
            y = y*y;
            y = z + y;
            for(int j=0; j<xdim; ++j)
            {
                x = (j - siz_x);
                x = x*x;
                x = sqrt(x + y);
                DIRECT_MULTIDIM_ELEM(testmap, n) = cos(2*PI/wavelength*x);
                ++n;
            }
        }
    }

    FileName fn;
    fn = formatString("fringes.vol");
    Image<double> saveImg;
    saveImg() = testmap;
    saveImg.write(fn);

    return testmap;
}

findCliffValue()

void Monogenic::findCliffValue	(	MultidimArray< double > &	inputmap,
		int &	radius,
		int &	radiuslimit,
		MultidimArray< int > &	mask,
		double	rsmooth
	)

Definition at line 67 of file monogenic_signal.cpp.

{
    double criticalZ = icdf_gauss(0.95);
    radiuslimit = XSIZE(inputmap)/2;
    double last_mean=0;
    double last_std2=1e-38;

    for (int rad = radius; rad<radiuslimit; rad++)
    {
        double sum=0;
        double sum2=0;
        double N=0;
        int sup;
        int inf;
        inf = rad*rad;
        sup = (rad + 1)*(rad + 1);
        FOR_ALL_ELEMENTS_IN_ARRAY3D(inputmap)
        {
            double aux = k*k + i*i + j*j;
            if ( (aux<sup) && (aux>=inf) )
            {
                double value = A3D_ELEM(inputmap, k, i, j);;
                sum += value;
                sum2 += value*value;
                N++;
            }
        }
        double mean = sum/N;
        double std2 = sum2/N - mean*mean;

        if (std2/last_std2<0.01)
        {
            radiuslimit = rad - 1;
            break;
        }

        last_mean = mean, last_std2 = std2;
    }

    std::cout << "There is no noise beyond a radius of " << radiuslimit << " px " << std::endl;
    std::cout << "Regions with a radius greater than " << radiuslimit << " px will not be considered" << std::endl;

    double raux = (radiuslimit - rsmooth);
    if (raux<=radius)
    {
        std::cout << "Warning: the boxsize is very close to "
                "the protein size please provide a greater box" << std::endl;
    }

    raux *= raux;

    FOR_ALL_ELEMENTS_IN_ARRAY3D(mask)
    {
        double aux = k*k + i*i + j*j;
        if (aux>=raux)
            A3D_ELEM(mask, k, i, j) = -1;
    }

}

fourierFreqs_3D()

MultidimArray< double > Monogenic::fourierFreqs_3D	(	const MultidimArray< std::complex< double > > &	myfftV,
		const MultidimArray< double > &	inputVol,
		Matrix1D< double > &	freq_fourier_x,
		Matrix1D< double > &	freq_fourier_y,
		Matrix1D< double > &	freq_fourier_z
	)

Definition at line 150 of file monogenic_signal.cpp.

{
    freq_fourier_z = fourierFreqVector(ZSIZE(myfftV), ZSIZE(inputVol));
        freq_fourier_y = fourierFreqVector(YSIZE(myfftV), YSIZE(inputVol));
        freq_fourier_x = fourierFreqVector(XSIZE(myfftV), XSIZE(inputVol));

    MultidimArray<double> iu;

    iu.initZeros(myfftV);

    double uz;
    double uy;
    double ux;
    double uz2;
    double u2;
    double uz2y2;
    long n=0;
    //  TODO: Take ZSIZE(myfftV) out of the loop
    //  TODO: Use freq_fourier_x instead of calling FFT_IDX2DIGFREQ

    for(size_t k=0; k<ZSIZE(myfftV); ++k)
    {
        FFT_IDX2DIGFREQ(k,ZSIZE(inputVol),uz);
        uz2 = uz*uz;
        for(size_t i=0; i<YSIZE(myfftV); ++i)
        {
            FFT_IDX2DIGFREQ(i,YSIZE(inputVol),uy);
            uz2y2 = uz2 + uy*uy;

            for(size_t j=0; j<XSIZE(myfftV); ++j)
            {
                FFT_IDX2DIGFREQ(j,XSIZE(inputVol), ux);
                u2 = uz2y2 + ux*ux;

                if ((k != 0) || (i != 0) || (j != 0))
                {
                    DIRECT_MULTIDIM_ELEM(iu,n) = 1/sqrt(u2);
                }
                else
                {
                    DIRECT_MULTIDIM_ELEM(iu,n) = 1e38;
                }
                ++n;
            }
        }
    }

    return iu;
}

fourierFreqVector()

Matrix1D< double > Monogenic::fourierFreqVector	(	size_t	dimarrayFourier,
		size_t	dimarrayReal
	)

Definition at line 132 of file monogenic_signal.cpp.

{
        double u;
        Matrix1D<double> freq_fourier;
          freq_fourier.initZeros(dimarrayFourier);
        VEC_ELEM(freq_fourier,0) = 1e-38; //A really low value to represent 0 avooiding singularities
    for(size_t k=1; k<dimarrayFourier; ++k){
        FFT_IDX2DIGFREQ(k,dimarrayReal, u);
        VEC_ELEM(freq_fourier, k) = u;
    }
    return freq_fourier;
}

monogenicAmplitude_3D_Fourier()

void Monogenic::monogenicAmplitude_3D_Fourier	(	const MultidimArray< std::complex< double > > &	myfftV,
		MultidimArray< double > &	iu,
		MultidimArray< double > &	amplitude,
		int	numberOfThreads
	)

Definition at line 575 of file monogenic_signal.cpp.

{
    Matrix1D<double> freq_fourier_z;
    Matrix1D<double> freq_fourier_y;
    Matrix1D<double> freq_fourier_x;

    iu = fourierFreqs_3D(myfftV, amplitude, freq_fourier_x, freq_fourier_y, freq_fourier_z);

    // Filter the input volume and add it to amplitude
    MultidimArray< std::complex<double> > fftVRiesz;
    MultidimArray< std::complex<double> > fftVRiesz_aux;
    fftVRiesz.initZeros(myfftV);
    fftVRiesz_aux.initZeros(myfftV);
    std::complex<double> J(0,1);

    long n=0;
    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(myfftV)
    {
        //double H=0.5*(1+cos((un-w1)*ideltal));
        DIRECT_MULTIDIM_ELEM(fftVRiesz, n) = DIRECT_MULTIDIM_ELEM(myfftV, n);
        DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n) = -J;
        DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n) *= DIRECT_MULTIDIM_ELEM(fftVRiesz, n);
        DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n) *= DIRECT_MULTIDIM_ELEM(iu, n);

    }
    MultidimArray<double> VRiesz;
    VRiesz.resizeNoCopy(amplitude);

    FourierTransformer transformer_inv;
        transformer_inv.setThreadsNumber(numberOfThreads);
    transformer_inv.inverseFourierTransform(fftVRiesz, VRiesz);

    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(amplitude)
        DIRECT_MULTIDIM_ELEM(amplitude,n) = DIRECT_MULTIDIM_ELEM(VRiesz,n)*DIRECT_MULTIDIM_ELEM(VRiesz,n);

    // Calculate first component of Riesz vector
    double uz;
    double uy;
    double ux;
    n=0;
    for(size_t k=0; k<ZSIZE(myfftV); ++k)
    {
        for(size_t i=0; i<YSIZE(myfftV); ++i)
        {
            for(size_t j=0; j<XSIZE(myfftV); ++j)
            {
                ux = VEC_ELEM(freq_fourier_x,j);
                DIRECT_MULTIDIM_ELEM(fftVRiesz, n) = ux*DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n);
                ++n;
            }
        }
    }

    transformer_inv.inverseFourierTransform(fftVRiesz, VRiesz);
    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(amplitude)
        DIRECT_MULTIDIM_ELEM(amplitude,n) += DIRECT_MULTIDIM_ELEM(VRiesz,n)*DIRECT_MULTIDIM_ELEM(VRiesz,n);

    // Calculate second and third components of Riesz vector
    n=0;
    for(size_t k=0; k<ZSIZE(myfftV); ++k)
    {
        uz = VEC_ELEM(freq_fourier_z,k);
        for(size_t i=0; i<YSIZE(myfftV); ++i)
        {
            uy = VEC_ELEM(freq_fourier_y,i);
            for(size_t j=0; j<XSIZE(myfftV); ++j)
            {
                DIRECT_MULTIDIM_ELEM(fftVRiesz, n) = uy*DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n);
                DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n) = uz*DIRECT_MULTIDIM_ELEM(fftVRiesz_aux, n);
                ++n;
            }
        }
    }
    transformer_inv.inverseFourierTransform(fftVRiesz, VRiesz);

    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(amplitude)
        DIRECT_MULTIDIM_ELEM(amplitude,n)+=DIRECT_MULTIDIM_ELEM(VRiesz,n)*DIRECT_MULTIDIM_ELEM(VRiesz,n);

    transformer_inv.inverseFourierTransform(fftVRiesz_aux, VRiesz);

    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(amplitude)
    {
        DIRECT_MULTIDIM_ELEM(amplitude,n) += DIRECT_MULTIDIM_ELEM(VRiesz,n)*DIRECT_MULTIDIM_ELEM(VRiesz,n);
        DIRECT_MULTIDIM_ELEM(amplitude,n) = sqrt(DIRECT_MULTIDIM_ELEM(amplitude,n));
    }
}

proteinRadiusVolumeAndShellStatistics()

void Monogenic::proteinRadiusVolumeAndShellStatistics	(	const MultidimArray< int > &	mask,
		int &	radius,
		long &	vol
	)

Definition at line 36 of file monogenic_signal.cpp.

{
    vol = 0;
    radius = 0;
    FOR_ALL_ELEMENTS_IN_ARRAY3D(mask)
    {

        if (A3D_ELEM(mask, k, i, j) == 1)
        {       
                        int R2 = (k*k + i*i + j*j);
            ++vol;
            if (R2>radius)
                radius = R2;
        }
    }
    radius = sqrt(radius);
    std::cout << "                                     " << std::endl;
    std::cout << "The protein has a radius of "<< radius << " px " << std::endl;
}

resolution2eval()

void Monogenic::resolution2eval	(	int &	count_res,
		double	step,
		double &	resolution,
		double &	last_resolution,
		double &	freq,
		double &	freqL,
		int &	last_fourier_idx,
		int &	volsize,
		bool &	continueIter,
		bool &	breakIter,
		double &	sampling,
		double &	maxRes
	)

Definition at line 212 of file monogenic_signal.cpp.

{
//TODO: simplify this function
    resolution = maxRes - count_res*step;

    freq = sampling/resolution;
    ++count_res;

    double Nyquist = 2*sampling;
    double aux_frequency;
    int fourier_idx;

    DIGFREQ2FFT_IDX(freq, volsize, fourier_idx);

    FFT_IDX2DIGFREQ(fourier_idx, volsize, aux_frequency);

    freq = aux_frequency;

    if (fourier_idx == last_fourier_idx)
    {
        continueIter = true;
        return;
    }

    last_fourier_idx = fourier_idx;
    resolution = sampling/aux_frequency;


    if (count_res == 0){
        last_resolution = resolution;
        }

    if (resolution<Nyquist)// || (resolution > last_resolution) )
    {
        breakIter = true;
        return;
    }

    freqL = sampling/(resolution + step);

    int fourier_idx_2;

    DIGFREQ2FFT_IDX(freqL, volsize, fourier_idx_2);

    if (fourier_idx_2 == fourier_idx)
    {
        if (fourier_idx > 0){
            FFT_IDX2DIGFREQ(fourier_idx - 1, volsize, freqL);
        }
        else{
            freqL = sampling/(resolution + step);
        }
    }

}

setLocalResolutionHalfMaps()

void Monogenic::setLocalResolutionHalfMaps	(	const MultidimArray< double > &	amplitudeMS,
		MultidimArray< int > &	pMask,
		MultidimArray< double > &	plocalResolutionMap,
		double	thresholdNoise,
		double	resolution,
		double	resolution_2
	)

Definition at line 518 of file monogenic_signal.cpp.

    {
        FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(amplitudeMS)
        {
            if (DIRECT_MULTIDIM_ELEM(pMask, n)>=1)
                if (DIRECT_MULTIDIM_ELEM(amplitudeMS, n)>thresholdNoise)
                {
                    DIRECT_MULTIDIM_ELEM(pMask, n) = 1;
                    DIRECT_MULTIDIM_ELEM(plocalResolutionMap, n) = resolution;
                }
                else{
                    DIRECT_MULTIDIM_ELEM(pMask, n) += 1;
                    if (DIRECT_MULTIDIM_ELEM(pMask, n) >2)
                    {
                        DIRECT_MULTIDIM_ELEM(pMask, n) = 0;
                        DIRECT_MULTIDIM_ELEM(plocalResolutionMap, n) = resolution_2;
                    }
                }
        }
    }

setLocalResolutionMap()

void Monogenic::setLocalResolutionMap	(	const MultidimArray< double > &	amplitudeMS,
		MultidimArray< int > &	pMask,
		MultidimArray< double > &	plocalResolutionMap,
		double	thresholdNoise,
		double	resolution,
		double	resolution_2
	)

Definition at line 548 of file monogenic_signal.cpp.

    {
        FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(amplitudeMS)
        {
            if (DIRECT_MULTIDIM_ELEM(pMask, n)>=1)
                if (DIRECT_MULTIDIM_ELEM(amplitudeMS, n)>thresholdNoise)
                {
                    DIRECT_MULTIDIM_ELEM(pMask, n) = 1;
                    DIRECT_MULTIDIM_ELEM(plocalResolutionMap, n) = resolution;//sampling/freq;
                }
                else{
                    DIRECT_MULTIDIM_ELEM(pMask, n) += 1;
                    if (DIRECT_MULTIDIM_ELEM(pMask, n) >2)
                    {
                        DIRECT_MULTIDIM_ELEM(pMask, n) = -1;
                        DIRECT_MULTIDIM_ELEM(plocalResolutionMap, n) = resolution_2;//maxRes - counter*R_;
                    }
                }
        }
    }

statisticsInBinaryMask2()

void Monogenic::statisticsInBinaryMask2	(	const MultidimArray< double > &	volS,
		const MultidimArray< double > &	volN,
		MultidimArray< int > &	mask,
		double &	meanS,
		double &	sdS2,
		double &	meanN,
		double &	sdN2,
		double &	significance,
		double &	thr95,
		double &	NS,
		double &	NN
	)

Definition at line 424 of file monogenic_signal.cpp.

{
    double sumS = 0;
    double sumS2 = 0;
    double sumN = 0;
    double sumN2 = 0;
    NN = 0;
    NS = 0;
    std::vector<double> noiseValues;

    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(volS)
    {
        if (DIRECT_MULTIDIM_ELEM(mask, n)>0)
        {
            double amplitudeValue=DIRECT_MULTIDIM_ELEM(volS, n);
            sumS  += amplitudeValue;
            sumS2 += amplitudeValue*amplitudeValue;
            ++NS;
        }
        if (DIRECT_MULTIDIM_ELEM(mask, n)>=0) //BE CAREFULL WITH THE =
        {
            double amplitudeValueN=DIRECT_MULTIDIM_ELEM(volN, n);
            noiseValues.push_back(amplitudeValueN);
            sumN  += amplitudeValueN;
            sumN2 += amplitudeValueN*amplitudeValueN;
            ++NN;
        }
    }

    std::sort(noiseValues.begin(),noiseValues.end());
    thr95 = noiseValues[size_t(noiseValues.size()*significance)];
    meanS = sumS/NS;
    meanN = sumN/NN;
    sdS2 = sumS2/NS - meanS*meanS;
    sdN2 = sumN2/NN - meanN*meanN;
}

statisticsInOutBinaryMask2()

void Monogenic::statisticsInOutBinaryMask2	(	const MultidimArray< double > &	volS,
		MultidimArray< int > &	mask,
		double &	meanS,
		double &	sdS2,
		double &	meanN,
		double &	sdN2,
		double &	significance,
		double &	thr95,
		double &	NS,
		double &	NN
	)

Definition at line 470 of file monogenic_signal.cpp.

{
    double sumS = 0;
    double sumS2 = 0;
    double sumN = 0;
    double sumN2 = 0;
    NN = 0;
    NS = 0;

    std::vector<double> noiseValues;

    FOR_ALL_DIRECT_ELEMENTS_IN_MULTIDIMARRAY(volS)
    {
        if (DIRECT_MULTIDIM_ELEM(mask, n)>=1)
        {
            double amplitudeValue=DIRECT_MULTIDIM_ELEM(volS, n);
            sumS  += amplitudeValue;
            sumS2 += amplitudeValue*amplitudeValue;
            ++NS;
        }
        if (DIRECT_MULTIDIM_ELEM(mask, n)==0)
        {
            double amplitudeValueN=DIRECT_MULTIDIM_ELEM(volS, n);
            noiseValues.push_back(amplitudeValueN);
            sumN  += amplitudeValueN;
            sumN2 += amplitudeValueN*amplitudeValueN;
            ++NN;
        }
    }

    std::sort(noiseValues.begin(),noiseValues.end());
    thr95 = noiseValues[size_t(noiseValues.size()*significance)];
    meanS = sumS/NS;
    meanN = sumN/NN;
    sdS2 = sumS2/NS - meanS*meanS;
    sdN2 = sumN2/NN - meanN*meanN;
}

---

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

• xmipp/libraries/data/monogenic_signal.h
• xmipp/libraries/data/monogenic_signal.cpp