#include <averager.h>
Inheritance diagram for EMAN::CtfAverager:
Public Member Functions | |
CtfAverager () | |
void | add_image (EMData *image) |
To add an image to the Averager. | |
EMData * | finish () |
Finish up the averaging and return the result. | |
vector< float > | get_snr () const |
Protected Attributes | |
XYData * | sf |
EMData * | curves |
bool | need_snr |
const char * | outfile |
Private Attributes | |
vector< float > | snr |
EMData * | image0_fft |
EMData * | image0_copy |
vector< vector< float > > | ctf |
vector< vector< float > > | ctfn |
float * | snri |
float * | snrn |
float * | tdr |
float * | tdi |
float * | tn |
float | filter |
int | nimg |
int | nx |
int | ny |
int | nz |
Definition at line 288 of file averager.h.
|
Definition at line 685 of file averager.cpp. 00685 : 00686 sf(0), curves(0), need_snr(false), outfile(0), 00687 image0_fft(0), image0_copy(0), snri(0), snrn(0), 00688 tdr(0), tdi(0), tn(0), 00689 filter(0), nimg(0), nx(0), ny(0), nz(0) 00690 { 00691 00692 }
|
|
To add an image to the Averager. This image will be averaged in this function.
Implements EMAN::Averager. Definition at line 694 of file averager.cpp. References EMAN::EMData::ap2ri(), EMAN::Ctf::apix, EMAN::Ctf::compute_1d(), EMAN::EMData::copy_head(), ctf, ctfn, curves, EMAN::EMData::do_fft(), EMAN::EMObject::f, EMAN::Util::fast_floor(), filter, EMAN::Dict::get(), EMAN::EMData::get_attr_dict(), EMAN::EMData::get_ctf(), EMAN::EMData::get_data(), EMAN::Averager::get_name(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), EMAN::EMData::get_zsize(), EMAN::EMData::has_ctff(), image0_copy, image0_fft, EMAN::EMUtil::is_same_size(), LOGERR, LOGWARN, nimg, nx, ny, EMAN::EMData::set_size(), snri, snrn, sqrt(), tdi, tdr, tn, EMAN::EMData::to_zero(), EMAN::EMData::update(), x, and y. 00695 { 00696 if (!image) { 00697 return; 00698 } 00699 00700 if (nimg >= 1 && !EMUtil::is_same_size(image, result)) { 00701 LOGERR("%sAverager can only process same-size Image", 00702 get_name().c_str()); 00703 return; 00704 } 00705 00706 if (image->get_zsize() != 1) { 00707 LOGERR("%sAverager: Only 2D images are currently supported", 00708 get_name().c_str()); 00709 } 00710 00711 string alg_name = get_name(); 00712 00713 if (alg_name == "CtfCW" || alg_name == "CtfCWauto") { 00714 if (image->get_ctf() != 0 && !image->has_ctff()) { 00715 LOGERR("%sAverager: Attempted CTF Correction with no ctf parameters", 00716 get_name().c_str()); 00717 } 00718 } 00719 else { 00720 if (image->get_ctf() != 0) { 00721 LOGERR("%sAverager: Attempted CTF Correction with no ctf parameters", 00722 get_name().c_str()); 00723 } 00724 } 00725 00726 nimg++; 00727 00728 00729 if (nimg == 1) { 00730 image0_fft = image->do_fft(); 00731 00732 nx = image0_fft->get_xsize(); 00733 ny = image0_fft->get_ysize(); 00734 nz = image0_fft->get_zsize(); 00735 00736 result = new EMData(); 00737 result->set_size(nx - 2, ny, nz); 00738 00739 00740 if (alg_name == "Weighting" && curves) { 00741 if (!sf) { 00742 LOGWARN("CTF curve in file will contain relative, not absolute SNR!"); 00743 } 00744 curves->set_size(Ctf::CTFOS * ny / 2, 3, 1); 00745 curves->to_zero(); 00746 } 00747 00748 00749 if (alg_name == "CtfC") { 00750 filter = params["filter"]; 00751 if (filter == 0) { 00752 filter = 22.0f; 00753 } 00754 float apix_y = image->get_attr_dict().get("apix_y"); 00755 float ds = 1.0f / (apix_y * ny * Ctf::CTFOS); 00756 filter = 1.0f / (filter * ds); 00757 } 00758 00759 if (alg_name == "CtfCWauto") { 00760 int nxy2 = nx * ny/2; 00761 00762 snri = new float[ny / 2]; 00763 snrn = new float[ny / 2]; 00764 tdr = new float[nxy2]; 00765 tdi = new float[nxy2]; 00766 tn = new float[nxy2]; 00767 00768 for (int i = 0; i < ny / 2; i++) { 00769 snri[i] = 0; 00770 snrn[i] = 0; 00771 } 00772 00773 for (int i = 0; i < nxy2; i++) { 00774 tdr[i] = 1; 00775 tdi[i] = 1; 00776 tn[i] = 1; 00777 } 00778 } 00779 00780 image0_copy = image0_fft->copy_head(); 00781 image0_copy->ap2ri(); 00782 image0_copy->to_zero(); 00783 } 00784 00785 Ctf::CtfType curve_type = Ctf::CTF_AMP; 00786 if (alg_name == "CtfCWauto") { 00787 curve_type = Ctf::CTF_AMP; 00788 } 00789 00790 float *src = image->get_data(); 00791 image->ap2ri(); 00792 Ctf *image_ctf = image->get_ctf(); 00793 int ny2 = image->get_ysize(); 00794 00795 vector<float> ctf1 = image_ctf->compute_1d(ny2,1.0f/(image_ctf->apix*image->get_ysize()), curve_type); 00796 00797 if (ctf1.size() == 0) { 00798 LOGERR("Unexpected CTF correction problem"); 00799 } 00800 00801 ctf.push_back(ctf1); 00802 00803 vector<float> ctfn1; 00804 if (sf) { 00805 ctfn1 = image_ctf->compute_1d(ny2,1.0f/(image_ctf->apix*image->get_ysize()), Ctf::CTF_SNR, sf); 00806 } 00807 else { 00808 ctfn1 = image_ctf->compute_1d(ny2,1.0f/(image_ctf->apix*image->get_ysize()), Ctf::CTF_SNR); 00809 } 00810 00811 ctfn.push_back(ctfn1); 00812 00813 if (alg_name == "CtfCWauto") { 00814 int j = 0; 00815 for (int y = 0; y < ny; y++) { 00816 for (int x = 0; x < nx / 2; x++, j += 2) { 00817 #ifdef _WIN32 00818 float r = (float)_hypot((float)x, (float)(y - ny / 2.0f)); 00819 #else 00820 float r = (float)hypot((float)x, (float)(y - ny / 2.0f)); 00821 #endif //_WIN32 00822 int l = static_cast < int >(Util::fast_floor(r)); 00823 00824 if (l >= 0 && l < ny / 2) { 00825 int k = y*nx/2 + x; 00826 tdr[k] *= src[j]; 00827 tdi[k] *= src[j + 1]; 00828 #ifdef _WIN32 00829 tn[k] *= (float)_hypot(src[j], src[j + 1]); 00830 #else 00831 tn[k] *= (float)hypot(src[j], src[j + 1]); 00832 #endif //_WIN32 00833 } 00834 } 00835 } 00836 } 00837 00838 00839 float *tmp_data = image0_copy->get_data(); 00840 00841 int j = 0; 00842 for (int y = 0; y < ny; y++) { 00843 for (int x = 0; x < nx / 2; x++, j += 2) { 00844 float r = Ctf::CTFOS * sqrt(x * x + (y - ny / 2.0f) * (y - ny / 2.0f)); 00845 int l = static_cast < int >(Util::fast_floor(r)); 00846 r -= l; 00847 00848 float f = 0; 00849 if (l <= Ctf::CTFOS * ny / 2 - 2) { 00850 f = (ctf1[l] * (1 - r) + ctf1[l + 1] * r); 00851 } 00852 tmp_data[j] += src[j] * f; 00853 tmp_data[j + 1] += src[j + 1] * f; 00854 } 00855 } 00856 00857 EMData *image_fft = image->do_fft(); 00858 image_fft->update(); 00859 if(image_ctf) {delete image_ctf; image_ctf=0;} 00860 }
|
|
Finish up the averaging and return the result.
Implements EMAN::Averager. Definition at line 862 of file averager.cpp. References ctf, ctfn, curves, EMAN::EMData::do_ift(), EMAN::Util::fast_floor(), filter, EMAN::EMData::get_data(), EMAN::Averager::get_name(), image0_copy, nimg, nx, ny, outfile, EMAN::Util::save_data(), snr, snri, snrn, tdi, tdr, tn, EMAN::EMData::update(), x, and y. 00863 { 00864 int j = 0; 00865 for (int y = 0; y < ny; y++) { 00866 for (int x = 0; x < nx / 2; x++, j += 2) { 00867 #ifdef _WIN32 00868 float r = (float) _hypot(x, y - ny / 2.0f); 00869 #else 00870 float r = (float) hypot(x, y - ny / 2.0f); 00871 #endif 00872 int l = static_cast < int >(Util::fast_floor(r)); 00873 if (l >= 0 && l < ny / 2) { 00874 int k = y*nx/2 + x; 00875 snri[l] += (tdr[k] + tdi[k]/tn[k]); 00876 snrn[l] += 1; 00877 } 00878 } 00879 } 00880 00881 for (int i = 0; i < ny / 2; i++) { 00882 snri[i] *= nimg / snrn[i]; 00883 } 00884 00885 if(strcmp(outfile, "") != 0) { 00886 Util::save_data(0, 1, snri, ny / 2, outfile); 00887 } 00888 00889 00890 float *cd = 0; 00891 if (curves) { 00892 cd = curves->get_data(); 00893 } 00894 00895 for (int i = 0; i < Ctf::CTFOS * ny / 2; i++) { 00896 float ctf0 = 0; 00897 for (int j = 0; j < nimg; j++) { 00898 ctf0 += ctfn[j][i]; 00899 if (ctf[j][i] == 0) { 00900 ctf[j][i] = 1.0e-12f; 00901 } 00902 00903 if (curves) { 00904 cd[i] += ctf[j][i] * ctfn[j][i]; 00905 cd[i + Ctf::CTFOS * ny / 2] += ctfn[j][i]; 00906 cd[i + 2 * Ctf::CTFOS * ny / 2] += ctfn[j][i]; 00907 } 00908 } 00909 00910 string alg_name = get_name(); 00911 00912 if (alg_name == "CtfCW" && need_snr) { 00913 snr[i] = ctf0; 00914 } 00915 00916 float ctf1 = ctf0; 00917 if (alg_name == "CtfCWauto") { 00918 ctf1 = snri[i / Ctf::CTFOS]; 00919 } 00920 00921 if (ctf1 <= 0.0001f) { 00922 ctf1 = 0.1f; 00923 } 00924 00925 if (alg_name == "CtfC") { 00926 for (int j = 0; j < nimg; j++) { 00927 ctf[j][i] = exp(-i * i / (filter * filter)) * ctfn[j][i] / (fabs(ctf[j][i]) * ctf1); 00928 } 00929 } 00930 else if (alg_name == "Weighting") { 00931 for (int j = 0; j < nimg; j++) { 00932 ctf[j][i] = ctfn[j][i] / ctf1; 00933 } 00934 } 00935 else if (alg_name == "CtfCW") { 00936 for (int j = 0; j < nimg; j++) { 00937 ctf[j][i] = (ctf1 / (ctf1 + 1)) * ctfn[j][i] / (ctf[j][i] * ctf1); 00938 } 00939 } 00940 else if (alg_name == "CtfCWauto") { 00941 for (int j = 0; j < nimg; j++) { 00942 ctf[j][i] = ctf1 * ctfn[j][i] / (fabs(ctf[j][i]) * ctf0); 00943 } 00944 } 00945 } 00946 00947 00948 if (curves) { 00949 for (int i = 0; i < Ctf::CTFOS * ny / 2; i++) { 00950 cd[i] /= cd[i + Ctf::CTFOS * ny / 2]; 00951 } 00952 curves->update(); 00953 } 00954 00955 image0_copy->update(); 00956 00957 float *result_data = result->get_data(); 00958 EMData *tmp_ift = image0_copy->do_ift(); 00959 float *tmp_ift_data = tmp_ift->get_data(); 00960 memcpy(result_data, tmp_ift_data, (nx - 2) * ny * sizeof(float)); 00961 00962 tmp_ift->update(); 00963 result->update(); 00964 00965 if( image0_copy ) 00966 { 00967 delete image0_copy; 00968 image0_copy = 0; 00969 } 00970 00971 if (snri) { 00972 delete[]snri; 00973 snri = 0; 00974 } 00975 00976 if (snrn) { 00977 delete[]snrn; 00978 snrn = 0; 00979 } 00980 00981 if( snri ) 00982 { 00983 delete [] snri; 00984 snri = 0; 00985 } 00986 if( snrn ) 00987 { 00988 delete [] snrn; 00989 snrn = 0; 00990 } 00991 if( tdr ) 00992 { 00993 delete [] tdr; 00994 tdr = 0; 00995 } 00996 if( tdi ) 00997 { 00998 delete [] tdi; 00999 tdi = 0; 01000 } 01001 if( tn ) 01002 { 01003 delete [] tn; 01004 tn = 0; 01005 } 01006 01007 return result; 01008 }
|
|
Definition at line 295 of file averager.h. 00296 {
00297 return snr;
00298 }
|
|
Definition at line 310 of file averager.h. Referenced by add_image(), and finish(). |
|
Definition at line 311 of file averager.h. Referenced by add_image(), and finish(). |
|
Definition at line 302 of file averager.h. Referenced by add_image(), and finish(). |
|
Definition at line 319 of file averager.h. Referenced by add_image(), and finish(). |
|
Definition at line 308 of file averager.h. Referenced by add_image(), and finish(). |
|
Definition at line 307 of file averager.h. Referenced by add_image(). |
|
Definition at line 303 of file averager.h. |
|
Definition at line 320 of file averager.h. Referenced by add_image(), and finish(). |
|
Definition at line 321 of file averager.h. |
|
Definition at line 322 of file averager.h. |
|
Definition at line 323 of file averager.h. |
|
Definition at line 304 of file averager.h. Referenced by finish(). |
|
Definition at line 301 of file averager.h. |
|
Definition at line 306 of file averager.h. Referenced by finish(). |
|
Definition at line 313 of file averager.h. Referenced by add_image(), and finish(). |
|
Definition at line 314 of file averager.h. Referenced by add_image(), and finish(). |
|
Definition at line 316 of file averager.h. Referenced by add_image(), and finish(). |
|
Definition at line 315 of file averager.h. Referenced by add_image(), and finish(). |
|
Definition at line 317 of file averager.h. Referenced by add_image(), and finish(). |