Both images should be FFTs. More...

#include <cmp.h>

Inheritance diagram for EMAN::TomoWedgeFscCmp:

Collaboration diagram for EMAN::TomoWedgeFscCmp:

[legend]

Public Member Functions
virtual float	cmp (EMData image, EMData with) const
	To compare 'image' with another image passed in through its parameters. More...

virtual string	get_name () const
	Get the Cmp's name. More...

virtual string	get_desc () const

TypeDict	get_param_types () const
	Get Cmp parameter information in a dictionary. More...

Public Member Functions inherited from EMAN::Cmp
virtual	~Cmp ()

virtual Dict	get_params () const
	Get the Cmp parameters in a key/value dictionary. More...

virtual void	set_params (const Dict &new_params)
	Set the Cmp parameters using a key/value dictionary. More...

Static Public Member Functions
static Cmp *	NEW ()

Static Public Attributes
static const string	NAME = "fsc.tomo.auto"

Additional Inherited Members
Protected Member Functions inherited from EMAN::Cmp
void	validate_input_args (const EMData image, const EMData with) const

Protected Attributes inherited from EMAN::Cmp
Dict	params

Detailed Description

Both images should be FFTs.

Computes a CCC excluding pixels with very small value (assumes they are missing wedge) Returns overlap amount as "fft_overlap" in "image"

Date: 2010-10-18

Definition at line 393 of file cmp.h.

Member Function Documentation

◆ cmp()

float TomoWedgeFscCmp::cmp	(	EMData *	image,
		EMData *	with
	)		const

virtual

To compare 'image' with another image passed in through its parameters.

An optional transformation may be used to transform the 2 images.

Parameters

image	The first image to be compared.
with	The second image to be comppared.

Returns: The comparison result. Smaller better by default

Implements EMAN::Cmp.

Definition at line 748 of file cmp.cpp.

{
        ENTERFUNC;
        EMData *imageo=image;
        EMData *witho=with;
        if (!image->is_complex()) image=image->do_fft();
        if (!with ->is_complex()) with =with ->do_fft();
        if (image->get_xsize()!=with->get_xsize() || image->get_ysize()!=with->get_ysize() || image->get_zsize()!=with->get_zsize())  throw InvalidCallException("Error: TomoWedgeFscCmp requires 2 images the same size");
 
        int nx=image->get_xsize();
        int ny=image->get_ysize();
        int nz=image->get_zsize();
 
        float sigmaimgval = params.set_default("sigmaimgval",0.5f);
        float sigmawithval = params.set_default("sigmawithval",0.5f);
        
        bool retcurve = params.set_default("retcurve",false);
        // User can pass in a sigma vector if they like, otherwise we call it 1/10 of the standard deviation in each shell
        // This has proven to be a good cutoff for identifying the missing wedge voxels, without throwing away too many
        // voxels if the image is complete in Fourier space (no wedge)
        vector<float> sigmaimg; 
        if (params.has_key("sigmaimg")) sigmaimg=params["sigmaimg"];
        else {
                sigmaimg=image->calc_radial_dist(nx/2,0,1,4);
                for (int i=0; i<nx/2; i++) sigmaimg[i]*=sigmaimg[i]*sigmaimgval;        // The value here is amplitude, we square to make comparison less expensive
        }
        
        vector<float> sigmawith;
        if (params.has_key("sigmawith")) sigmawith = params["sigmawith"];
        else {
                sigmawith=with->calc_radial_dist(nx/2,0,1,4);
                for (int i=0; i<nx/2; i++) sigmawith[i]*=sigmawith[i]*sigmawithval; // The value here is amplitude, we square to make comparison less expensive
        }
 
        float apix = params.set_default("apix",image->get_attr_default("apix_x", 1.0f));
        //get min and max res
        float minres = params.set_default("minres",0.0f);
        float maxres = params.set_default("maxres", 0.0f);
        
        int pmin = params.set_default("pmin",3);
        int pmax = params.set_default("pmax", ny/2);
        
        if (minres>0) 
                pmin=(int)floor((apix*ny)/minres);              //cutoff in pixels, assume square
        if (maxres>0) 
                pmax=(int)ceil((apix*ny)/maxres);
 
        int negative = params.set_default("negative",1);
        
        double sum=0;
        double sumsq1=0;
        double sumsq2=0;
        double norm=0;
        vector<float> curve(nx/2), sqcv1(nx/2), sqcv2(nx/2);
        for (int z=0; z<nz; z++) {
                int za=z<nz/2?z:nz-z;
                if (za>pmax) continue;
                for (int y=0; y<ny; y++) {
                        int ya=y<ny/2?y:ny-y;
                        if (ya>pmax) continue;
                        for (int x=0; x<nx; x+=2) {
                                float r2=Util::hypot3sq(x/2,ya,za);     // origin at 0,0; periodic
                                int r=int(sqrtf(r2));
//                              float rf=Util::hypot3(x/2,y<ny/2?y:ny-y,z<nz/2?z:nz-z); // origin at 0,0; periodic
//                              int r=int(rf);
                                if (r<pmin || r>pmax) continue;
                                
                                float v1r=image->get_value_at(x,y,z);
                                float v1i=image->get_value_at(x+1,y,z);
                                float v1=Util::square_sum(v1r,v1i);
                                if (v1<sigmaimg[r]) continue;
                                
                                float v2r=with->get_value_at(x,y,z);
                                float v2i=with->get_value_at(x+1,y,z);
                                float v2=Util::square_sum(v2r,v2i);
                                if (v2<sigmawith[r]) continue;
                                
                                sum+=(v1r*v2r+v1i*v2i)/r2;      // The r downweighting of points allows us to integrate over the entire image rather than computing a FSC and integrating that
                                sumsq1+=v1/r2;
//                              if (Util::is_nan(sumsq1)) { printf("TomoWedgeCccCmp: NaN encountered: %d %d %d %f %f %f\n",x,y,z,v1r,v1i,v1); }
                                sumsq2+=v2/r2;
                                norm+=1.0;
                                if (retcurve){
                                        curve[r]+=v1r*v2r+v1i*v2i;
                                        sqcv1[r]+=v1;
                                        sqcv2[r]+=v2;
                                }
                        }
                }
        }
        image->set_attr("fft_overlap",(float)(2.0*norm/(image->get_xsize()*image->get_ysize()*image->get_zsize())));
//      printf("%f\t%f\t%f\t%f\t%f\n",s1,s2,sumsq1,sumsq2,norm);
        if (retcurve){
                for(int i=0; i<nx/2; i++){
                        if (sqcv1[i]==0) sqcv1[i]=1;
                        if (sqcv2[i]==0) sqcv2[i]=1;
                        curve[i]=curve[i]/(sqrt(sqcv1[i])*sqrt(sqcv2[i]));
                }
                image->set_attr("fsc_curve", curve);
        }
 
        if (imageo!=image) delete image;
        if (witho!=with) delete with;
        
        if (negative) sum*=-1.0;
        return float(sum/(sqrt(sumsq1)*sqrt(sumsq2)));
}

References EMAN::EMData::calc_radial_dist(), ENTERFUNC, EMAN::Dict::has_key(), EMAN::Util::hypot3sq(), InvalidCallException, is_complex(), EMAN::Cmp::params, EMAN::Dict::set_default(), sqrt(), EMAN::Util::square_sum(), x, and y.

◆ get_desc()

virtual string EMAN::TomoWedgeFscCmp::get_desc ( ) const

inlinevirtual

Implements EMAN::Cmp.

Definition at line 403 of file cmp.h.

                {
                        return "FSC of two FFTs with missing wedge taken into account via a dynamic threshold. Also stores overlap as 'fft_overlap' in \
first input image on return. This comparator can take the missing wedge into account automatically without any precomputation of wedge \
orientation. If sigmaimg and/or sigmawith are provided, they must contain floating point arrays with at least as many elements as the \
Fourier radius in pixels, and these curves represent the thresholds for each image for including the values in the FSC, or considering them \
to be missing values. If not provided, these values will be computed automatically. These options are primarily provided for speed.";
                }

◆ get_name()

virtual string EMAN::TomoWedgeFscCmp::get_name ( ) const

inlinevirtual

Get the Cmp's name.

Each Cmp is identified by a unique name.

Returns: The Cmp's name.

Implements EMAN::Cmp.

Definition at line 398 of file cmp.h.

                {
                        return NAME;
                }

References NAME.

◆ get_param_types()

TypeDict EMAN::TomoWedgeFscCmp::get_param_types ( ) const

inlinevirtual

Get Cmp parameter information in a dictionary.

Each parameter has one record in the dictionary. Each record contains its name, data-type, and description.

Returns: A dictionary containing the parameter info.

Implements EMAN::Cmp.

Definition at line 417 of file cmp.h.

                {
                        TypeDict d;
//                      d.put("norm", EMObject::BOOL,"Whether the cross correlation image should be normalized (should be for normalized images). Default is true.");
//                      d.put("ccf", EMObject::EMDATA,"The ccf image, can be provided if it already exists to avoid recalculating it");
                        d.put("sigmaimg", EMObject::FLOATARRAY, "Res dependent amplitude threshold. If provided this must contain a floating point array with as many elements as the Fourier radius of the image. ");
                        d.put("sigmawith", EMObject::FLOATARRAY, "Resolution dependent coefficient for thresholding values included in the dot product in the 'with' image. Default = 0.1 and is normally fine");
                        d.put("sigmaimgval", EMObject::FLOAT, "Sigma coefficient for thresholding values included in the dot product. default = 0.5");
                        d.put("sigmawithval", EMObject::FLOAT, "Sigma coefficient for thresholding values included in the dot product in the 'with' image. Default = 0.5");
                        d.put("minres", EMObject::FLOAT, "The minimum resolution to accept (1/A) Default is 3 pixels. overwrites pmin");
                        d.put("maxres", EMObject::FLOAT, "The maximum resolution to accept (1/A) Default is axial Nyquist. overwrites pmax");
                        d.put("pmin", EMObject::INT, "The minimum resolution in pixels.");
                        d.put("pmax", EMObject::INT, "The maximum resolution in pixels.");
                        d.put("retcurve", EMObject::BOOL, "Return fsc curve in image header");
                        return d;
                }

References EMAN::EMObject::BOOL, EMAN::EMObject::FLOAT, EMAN::EMObject::FLOATARRAY, EMAN::EMObject::INT, and EMAN::TypeDict::put().

◆ NEW()

static Cmp * EMAN::TomoWedgeFscCmp::NEW ( )

inlinestatic

Definition at line 412 of file cmp.h.

                {
                        return new TomoWedgeFscCmp();
                }

Member Data Documentation

◆ NAME

const string TomoWedgeFscCmp::NAME = "fsc.tomo.auto"

static

Definition at line 434 of file cmp.h.

Referenced by get_name().

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

Public Member Functions

Static Public Member Functions

Static Public Attributes

Additional Inherited Members

Detailed Description

Member Function Documentation

◆ cmp()

◆ get_desc()

◆ get_name()

◆ get_param_types()

◆ NEW()

Member Data Documentation

◆ NAME