00001
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036 #include "averager.h"
00037 #include "emdata.h"
00038 #include "xydata.h"
00039 #include "ctf.h"
00040 #include <cstring>
00041 #include "plugins/averager_template.h"
00042
00043 using namespace EMAN;
00044
00045 const string ImageAverager::NAME = "mean";
00046 const string TomoAverager::NAME = "mean.tomo";
00047 const string MinMaxAverager::NAME = "minmax";
00048 const string AbsMaxMinAverager::NAME = "absmaxmin";
00049 const string IterationAverager::NAME = "iteration";
00050 const string WeightingAverager::NAME = "snr_weight";
00051 const string CtfCAverager::NAME = "ctfc";
00052 const string CtfCWAverager::NAME = "ctfcw";
00053 const string CtfCWautoAverager::NAME = "ctfw.auto";
00054 const string CtfCAutoAverager::NAME = "ctf.auto";
00055
00056 template <> Factory < Averager >::Factory()
00057 {
00058 force_add<ImageAverager>();
00059 force_add<MinMaxAverager>();
00060 force_add<AbsMaxMinAverager>();
00061 force_add<IterationAverager>();
00062
00063
00064
00065
00066 force_add<CtfCWautoAverager>();
00067 force_add<CtfCAutoAverager>();
00068 force_add<TomoAverager>();
00069
00070 }
00071
00072 void Averager::mult(const float& s)
00073 {
00074 if ( result != 0 )
00075 {
00076 result->mult(s);
00077 }
00078 else throw NullPointerException("Error, attempted to multiply the result image, which is NULL");
00079 }
00080
00081
00082 void Averager::add_image_list(const vector<EMData*> & image_list)
00083 {
00084 for (size_t i = 0; i < image_list.size(); i++) {
00085 add_image(image_list[i]);
00086 }
00087 }
00088
00089 TomoAverager::TomoAverager()
00090 : norm_image(0)
00091 {
00092
00093 }
00094
00095 void TomoAverager::add_image(EMData * image)
00096 {
00097 if (!image) {
00098 return;
00099 }
00100
00101 if (!image->is_complex()) {
00102 image=image->do_fft();
00103 image->set_attr("free_me",(int)1);
00104 }
00105
00106 if (result!=0 && !EMUtil::is_same_size(image, result)) {
00107 LOGERR("%s Averager can only process same-size Images",
00108 get_name().c_str());
00109 return;
00110 }
00111
00112 int nx = image->get_xsize();
00113 int ny = image->get_ysize();
00114 int nz = image->get_zsize();
00115 size_t image_size = (size_t)nx * ny * nz;
00116
00117 if (norm_image == 0) {
00118 printf("init average %d %d %d",nx,ny,nz);
00119 result = image->copy_head();
00120 result->to_zero();
00121
00122 norm_image = image->copy_head();
00123 norm_image->to_zero();
00124
00125 thresh_sigma = (float)params.set_default("thresh_sigma", 0.01);
00126 }
00127
00128 float *result_data = result->get_data();
00129 float *norm_data = norm_image->get_data();
00130 float *data = image->get_data();
00131 float thresh=image->get_attr("sigma");
00132 thresh=2.0f*thresh*thresh*thresh_sigma;
00133
00134
00135 for (size_t j = 0; j < image_size; j+=2) {
00136 float f=data[j];
00137 float g=data[j+1];
00138 float inten=f*f+g*g;
00139
00140 if (inten<thresh) continue;
00141
00142 result_data[j] +=f;
00143 result_data[j+1]+=g;
00144
00145 norm_data[j] +=1.0;
00146 norm_data[j+1]+=1.0;
00147 }
00148
00149 if (image->has_attr("free_me")) delete image;
00150 }
00151
00152 EMData * TomoAverager::finish()
00153 {
00154 if (norm_image==0 || result==0) return NULL;
00155
00156 int nx = result->get_xsize();
00157 int ny = result->get_ysize();
00158 int nz = result->get_zsize();
00159 size_t image_size = (size_t)nx * ny * nz;
00160
00161 float *result_data = result->get_data();
00162 float *norm_data = norm_image->get_data();
00163
00164 printf("finish average %d %d %d",nx,ny,nz);
00165
00166 for (size_t j = 0; j < image_size; j++) {
00167 if (norm_data[j]==0.0) result_data[j]=0.0;
00168 else result_data[j]/=norm_data[j];
00169 }
00170
00171 norm_image->update();
00172 result->update();
00173
00174 EMData *ret = result->do_ift();
00175 ret->set_attr("ptcl_repr",norm_image->get_attr("maximum"));
00176 if ((int)params.set_default("save_norm", 0))
00177 norm_image->write_image("norm.hdf");
00178
00179 delete result;
00180 delete norm_image;
00181 result=0;
00182 norm_image=0;
00183
00184 return ret;
00185 }
00186
00187
00188 ImageAverager::ImageAverager()
00189 : sigma_image(0), nimg_n0(0), ignore0(0), nimg(0)
00190 {
00191
00192 }
00193
00194 void ImageAverager::add_image(EMData * image)
00195 {
00196 if (!image) {
00197 return;
00198 }
00199
00200 if (nimg >= 1 && !EMUtil::is_same_size(image, result)) {
00201 LOGERR("%sAverager can only process same-size Image",
00202 get_name().c_str());
00203 return;
00204 }
00205
00206 nimg++;
00207
00208 int nx = image->get_xsize();
00209 int ny = image->get_ysize();
00210 int nz = image->get_zsize();
00211 size_t image_size = (size_t)nx * ny * nz;
00212
00213 if (nimg == 1) {
00214 result = image->copy_head();
00215 result->set_size(nx, ny, nz);
00216 sigma_image = params.set_default("sigma", (EMData*)0);
00217 ignore0 = params["ignore0"];
00218
00219 nimg_n0 = new int[image_size];
00220 for (size_t i = 0; i < image_size; ++i) {
00221 nimg_n0[i] = 0;
00222 }
00223 }
00224
00225 float *result_data = result->get_data();
00226 float *sigma_image_data = 0;
00227 if (sigma_image) {
00228 sigma_image->set_size(nx, ny, nz);
00229 sigma_image_data = sigma_image->get_data();
00230 }
00231
00232 float * image_data = image->get_data();
00233
00234 if (!ignore0) {
00235 for (size_t j = 0; j < image_size; ++j) {
00236 float f = image_data[j];
00237 result_data[j] += f;
00238 if (sigma_image_data) {
00239 sigma_image_data[j] += f * f;
00240 }
00241 }
00242 }
00243 else {
00244 for (size_t j = 0; j < image_size; ++j) {
00245 float f = image_data[j];
00246 if (f) {
00247 result_data[j] += f;
00248 if (sigma_image_data) {
00249 sigma_image_data[j] += f * f;
00250 }
00251 nimg_n0[j]++;
00252 }
00253 }
00254 }
00255 }
00256
00257 EMData * ImageAverager::finish()
00258 {
00259 if (result && nimg > 1) {
00260 size_t image_size = (size_t)result->get_xsize() * result->get_ysize() * result->get_zsize();
00261 float * result_data = result->get_data();
00262
00263 if (!ignore0) {
00264 for (size_t j = 0; j < image_size; ++j) {
00265 result_data[j] /= nimg;
00266 }
00267
00268 if (sigma_image) {
00269 float * sigma_image_data = sigma_image->get_data();
00270
00271 for (size_t j = 0; j < image_size; ++j) {
00272 float f1 = sigma_image_data[j] / nimg;
00273 float f2 = result_data[j];
00274 sigma_image_data[j] = sqrt(f1 - f2 * f2);
00275 }
00276
00277 sigma_image->update();
00278 }
00279 }
00280 else {
00281 for (size_t j = 0; j < image_size; ++j) {
00282 if (nimg_n0[j]>0) result_data[j] /= nimg_n0[j];
00283 }
00284 if (sigma_image) {
00285 float * sigma_image_data = sigma_image->get_data();
00286
00287 for (size_t j = 0; j < image_size; ++j) {
00288 float f1 = sigma_image_data[j] / nimg_n0[j];
00289 float f2 = result_data[j];
00290 sigma_image_data[j] = sqrt(f1 - f2 * f2);
00291 }
00292
00293 sigma_image->update();
00294 }
00295 }
00296
00297 result->update();
00298
00299 }
00300 result->set_attr("ptcl_repr",nimg);
00301
00302 if( nimg_n0 )
00303 {
00304 delete [] nimg_n0;
00305 nimg_n0 = 0;
00306 }
00307
00308 return result;
00309 }
00310
00311 #if 0
00312 EMData *ImageAverager::average(const vector < EMData * >&image_list) const
00313 {
00314 if (image_list.size() == 0) {
00315 return 0;
00316 }
00317
00318 EMData *sigma_image = params["sigma"];
00319 int ignore0 = params["ignore0"];
00320
00321 EMData *image0 = image_list[0];
00322
00323 int nx = image0->get_xsize();
00324 int ny = image0->get_ysize();
00325 int nz = image0->get_zsize();
00326 size_t image_size = (size_t)nx * ny * nz;
00327
00328 EMData *result = new EMData();
00329 result->set_size(nx, ny, nz);
00330
00331 float *result_data = result->get_data();
00332 float *sigma_image_data = 0;
00333
00334 if (sigma_image) {
00335 sigma_image->set_size(nx, ny, nz);
00336 sigma_image_data = sigma_image->get_data();
00337 }
00338
00339 int c = 1;
00340 if (ignore0) {
00341 for (size_t j = 0; j < image_size; ++j) {
00342 int g = 0;
00343 for (size_t i = 0; i < image_list.size(); i++) {
00344 float f = image_list[i]->get_value_at(j);
00345 if (f) {
00346 g++;
00347 result_data[j] += f;
00348 if (sigma_image_data) {
00349 sigma_image_data[j] += f * f;
00350 }
00351 }
00352 }
00353 if (g) {
00354 result_data[j] /= g;
00355 }
00356 }
00357 }
00358 else {
00359 float *image0_data = image0->get_data();
00360 if (sigma_image_data) {
00361 memcpy(sigma_image_data, image0_data, image_size * sizeof(float));
00362
00363 for (size_t j = 0; j < image_size; ++j) {
00364 sigma_image_data[j] *= sigma_image_data[j];
00365 }
00366 }
00367
00368 image0->update();
00369 memcpy(result_data, image0_data, image_size * sizeof(float));
00370
00371 for (size_t i = 1; i < image_list.size(); i++) {
00372 EMData *image = image_list[i];
00373
00374 if (EMUtil::is_same_size(image, result)) {
00375 float *image_data = image->get_data();
00376
00377 for (size_t j = 0; j < image_size; ++j) {
00378 result_data[j] += image_data[j];
00379 }
00380
00381 if (sigma_image_data) {
00382 for (size_t j = 0; j < image_size; ++j) {
00383 sigma_image_data[j] += image_data[j] * image_data[j];
00384 }
00385 }
00386
00387 image->update();
00388 c++;
00389 }
00390 }
00391
00392 for (size_t j = 0; j < image_size; ++j) {
00393 result_data[j] /= static_cast < float >(c);
00394 }
00395 }
00396
00397 if (sigma_image_data) {
00398 for (size_t j = 0; j < image_size; ++j) {
00399 float f1 = sigma_image_data[j] / c;
00400 float f2 = result_data[j];
00401 sigma_image_data[j] = sqrt(f1 - f2 * f2);
00402 }
00403 }
00404
00405 if (sigma_image_data) {
00406 sigma_image->update();
00407 }
00408
00409 result->update();
00410 return result;
00411 }
00412 #endif
00413
00414 MinMaxAverager::MinMaxAverager()
00415 : nimg(0)
00416 {
00417
00418
00419 max = 0;
00420 }
00421
00422 void MinMaxAverager::add_image(EMData * image)
00423 {
00424 if (!image) {
00425 return;
00426 }
00427
00428 if (nimg >= 1 && !EMUtil::is_same_size(image, result)) {
00429 LOGERR("%sAverager can only process same-size Image",
00430 get_name().c_str());
00431 return;
00432 }
00433
00434 nimg++;
00435
00436 int nx = image->get_xsize();
00437 int ny = image->get_ysize();
00438 int nz = image->get_zsize();
00439
00440 if (nimg == 1) {
00441 result = image->copy();
00442 max = params["max"];
00443 return;
00444 }
00445
00446 for (int z=0; z<nz; z++) {
00447 for (int y=0; y<ny; y++) {
00448 for (int x=0; x<nx; x++) {
00449 if (result->get_value_at(x,y,z)>image->get_value_at(x,y,z))
00450 { if (!max) result->set_value_at(x,y,z,image->get_value_at(x,y,z)); }
00451 else { if (max) result->set_value_at(x,y,z,image->get_value_at(x,y,z)); }
00452 }
00453 }
00454 }
00455
00456 }
00457
00458 EMData *MinMaxAverager::finish()
00459 {
00460 result->update();
00461 result->set_attr("ptcl_repr",nimg);
00462
00463 if (result && nimg > 1) return result;
00464
00465 return NULL;
00466 }
00467
00468 AbsMaxMinAverager::AbsMaxMinAverager() : nimg(0)
00469 {
00470
00471
00472 min = 0;
00473 }
00474
00475 void AbsMaxMinAverager::add_image(EMData * image)
00476 {
00477 if (!image) {
00478 return;
00479 }
00480
00481 if (nimg >= 1 && !EMUtil::is_same_size(image, result)) {
00482 LOGERR("%sAverager can only process same-size Image",
00483 get_name().c_str());
00484 return;
00485 }
00486
00487 nimg++;
00488
00489 int nx = image->get_xsize();
00490 int ny = image->get_ysize();
00491 int nz = image->get_zsize();
00492
00493 size_t imgsize = (size_t)nx*ny*nz;
00494
00495 if (nimg == 1) {
00496 result = image->copy();
00497 min = params["min"];
00498 return;
00499 }
00500
00501 float * data = result->get_data();
00502 float * src_data = image->get_data();
00503
00504 for(size_t i=0; i<imgsize; ++i) {
00505 if(!min) {
00506 if (fabs(data[i]) < fabs(src_data[i])) data[i]=src_data[i];
00507 }
00508 else {
00509 if (fabs(data[i]) > fabs(src_data[i])) data[i]=src_data[i];
00510 }
00511 }
00512 }
00513
00514 EMData *AbsMaxMinAverager::finish()
00515 {
00516 result->update();
00517 result->set_attr("ptcl_repr",nimg);
00518
00519 if (result && nimg > 1) return result;
00520
00521 return NULL;
00522 }
00523
00524 IterationAverager::IterationAverager() : nimg(0)
00525 {
00526
00527 }
00528
00529 void IterationAverager::add_image( EMData * image)
00530 {
00531 if (!image) {
00532 return;
00533 }
00534
00535 if (nimg >= 1 && !EMUtil::is_same_size(image, result)) {
00536 LOGERR("%sAverager can only process same-size Image",
00537 get_name().c_str());
00538 return;
00539 }
00540
00541 nimg++;
00542
00543 int nx = image->get_xsize();
00544 int ny = image->get_ysize();
00545 int nz = image->get_zsize();
00546 size_t image_size = (size_t)nx * ny * nz;
00547
00548 if (nimg == 1) {
00549 result = image->copy_head();
00550 result->set_size(nx, ny, nz);
00551 sigma_image = image->copy_head();
00552 sigma_image->set_size(nx, ny, nz);
00553 }
00554
00555 float *image_data = image->get_data();
00556 float *result_data = result->get_data();
00557 float *sigma_image_data = sigma_image->get_data();
00558
00559 for (size_t j = 0; j < image_size; ++j) {
00560 float f = image_data[j];
00561 result_data[j] += f;
00562 sigma_image_data[j] += f * f;
00563 }
00564
00565
00566 }
00567
00568 EMData * IterationAverager::finish()
00569 {
00570 if (nimg < 1) {
00571 return result;
00572 }
00573
00574 int nx = result->get_xsize();
00575 int ny = result->get_ysize();
00576 int nz = result->get_zsize();
00577 size_t image_size = (size_t)nx * ny * nz;
00578
00579 float *result_data = result->get_data();
00580 float *sigma_image_data = sigma_image->get_data();
00581
00582 for (size_t j = 0; j < image_size; ++j) {
00583 result_data[j] /= nimg;
00584 float f1 = sigma_image_data[j] / nimg;
00585 float f2 = result_data[j];
00586 sigma_image_data[j] = sqrt(f1 - f2 * f2) / sqrt((float)nimg);
00587 }
00588
00589 result->update();
00590 sigma_image->update();
00591
00592 result->append_image("iter.hed");
00593 float sigma = sigma_image->get_attr("sigma");
00594 float *sigma_image_data2 = sigma_image->get_data();
00595 float *result_data2 = result->get_data();
00596 float *d2 = new float[nx * ny];
00597 size_t sec_size = nx * ny * sizeof(float);
00598
00599 memcpy(d2, result_data2, sec_size);
00600 memcpy(sigma_image_data2, result_data2, sec_size);
00601
00602 printf("Iter sigma=%f\n", sigma);
00603
00604 for (int k = 0; k < 1000; k++) {
00605 for (int i = 1; i < nx - 1; i++) {
00606 for (int j = 1; j < ny - 1; j++) {
00607 int l = i + j * nx;
00608 float c1 = (d2[l - 1] + d2[l + 1] + d2[l - nx] + d2[l + nx]) / 4.0f - d2[l];
00609 float c2 = fabs(result_data2[l] - sigma_image_data2[l]) / sigma;
00610 result_data2[l] += c1 * Util::eman_erfc(c2) / 100.0f;
00611 }
00612 }
00613
00614 memcpy(d2, result_data2, sec_size);
00615 }
00616
00617 if( d2 )
00618 {
00619 delete[]d2;
00620 d2 = 0;
00621 }
00622
00623 sigma_image->update();
00624 if( sigma_image )
00625 {
00626 delete sigma_image;
00627 sigma_image = 0;
00628 }
00629
00630 result->update();
00631 result->append_image("iter.hed");
00632
00633
00634 return result;
00635 }
00636
00637 CtfCWautoAverager::CtfCWautoAverager()
00638 : nimg(0)
00639 {
00640
00641 }
00642
00643
00644 void CtfCWautoAverager::add_image(EMData * image)
00645 {
00646 if (!image) {
00647 return;
00648 }
00649
00650
00651
00652 EMData *fft=image->do_fft();
00653
00654 if (nimg >= 1 && !EMUtil::is_same_size(fft, result)) {
00655 LOGERR("%s Averager can only process images of the same size", get_name().c_str());
00656 return;
00657 }
00658
00659 nimg++;
00660 if (nimg == 1) {
00661 result = fft->copy_head();
00662 result->to_zero();
00663 }
00664
00665 Ctf *ctf = (Ctf *)image->get_attr("ctf");
00666
00667
00668
00669
00670 float b=ctf->bfactor;
00671 ctf->bfactor=100.0;
00672
00673
00674
00675 EMData *snr = result -> copy();
00676 ctf->compute_2d_complex(snr,Ctf::CTF_SNR);
00677
00678 EMData *ctfi = result-> copy();
00679 ctf->compute_2d_complex(ctfi,Ctf::CTF_AMP);
00680
00681 ctf->bfactor=b;
00682
00683 float *outd = result->get_data();
00684 float *ind = fft->get_data();
00685 float *snrd = snr->get_data();
00686 float *ctfd = ctfi->get_data();
00687
00688 size_t sz=snr->get_xsize()*snr->get_ysize();
00689 for (size_t i = 0; i < sz; i+=2) {
00690 if (snrd[i]<0) snrd[i]=0;
00691 ctfd[i]=fabs(ctfd[i]);
00692 if (ctfd[i]<.05) ctfd[i]=0.05f;
00693
00694
00695
00696
00697 outd[i]+=ind[i]*snrd[i]/ctfd[i];
00698 outd[i+1]+=ind[i+1]*snrd[i]/ctfd[i];
00699 }
00700
00701 if (nimg==1) {
00702 snrsum=snr->copy_head();
00703 float *ssnrd=snrsum->get_data();
00704
00705 for (size_t i = 0; i < sz; i+=2) { ssnrd[i]=1.0; ssnrd[i+1]=0.0; }
00706 }
00707
00708 snr->process_inplace("math.absvalue");
00709 snrsum->add(*snr);
00710
00711 delete ctf;
00712 delete fft;
00713 delete snr;
00714 delete ctfi;
00715 }
00716
00717 EMData * CtfCWautoAverager::finish()
00718 {
00719
00720
00721
00722
00723
00724
00725
00726
00727
00728
00729 int nx=result->get_xsize();
00730 int ny=result->get_ysize();
00731 float *snrsd=snrsum->get_data();
00732 float *outd=result->get_data();
00733
00734 int rm=(ny-2)*(ny-2)/4;
00735 for (int j=0; j<ny; j++) {
00736 for (int i=0; i<nx; i+=2) {
00737 size_t ii=i+j*nx;
00738 if ((j<ny/2 && i*i/4+j*j>rm) ||(j>=ny/2 && i*i/4+(ny-j)*(ny-j)>rm) || snrsd[ii]==0) { outd[ii]=outd[ii+1]=0; continue; }
00739 outd[ii]/=snrsd[ii];
00740 outd[ii+1]/=snrsd[ii];
00741 }
00742 }
00743
00744 result->update();
00745 result->set_attr("ptcl_repr",nimg);
00746 result->set_attr("ctf_snr_total",snrsum->calc_radial_dist(snrsum->get_ysize()/2,0,1,false));
00747 result->set_attr("ctf_wiener_filtered",1);
00748
00749 delete snrsum;
00750 EMData *ret=result->do_ift();
00751 delete result;
00752 result=NULL;
00753 return ret;
00754 }
00755
00756 CtfCAutoAverager::CtfCAutoAverager()
00757 : nimg(0)
00758 {
00759
00760 }
00761
00762
00763 void CtfCAutoAverager::add_image(EMData * image)
00764 {
00765 if (!image) {
00766 return;
00767 }
00768
00769
00770
00771 EMData *fft=image->do_fft();
00772
00773 if (nimg >= 1 && !EMUtil::is_same_size(fft, result)) {
00774 LOGERR("%s Averager can only process images of the same size", get_name().c_str());
00775 return;
00776 }
00777
00778 nimg++;
00779 if (nimg == 1) {
00780 result = fft->copy_head();
00781 result->to_zero();
00782 }
00783
00784 Ctf *ctf = (Ctf *)image->get_attr("ctf");
00785 float b=ctf->bfactor;
00786 ctf->bfactor=0;
00787
00788 EMData *snr = result -> copy();
00789 ctf->compute_2d_complex(snr,Ctf::CTF_SNR);
00790 EMData *ctfi = result-> copy();
00791 ctf->compute_2d_complex(ctfi,Ctf::CTF_AMP);
00792
00793 ctf->bfactor=b;
00794
00795 float *outd = result->get_data();
00796 float *ind = fft->get_data();
00797 float *snrd = snr->get_data();
00798 float *ctfd = ctfi->get_data();
00799
00800 size_t sz=snr->get_xsize()*snr->get_ysize();
00801 for (size_t i = 0; i < sz; i+=2) {
00802 if (snrd[i]<0) snrd[i]=0;
00803 ctfd[i]=fabs(ctfd[i]);
00804
00805
00806 if (ctfd[i]<.05) ctfd[i]=0.05f;
00807
00808
00809
00810
00811
00812
00813 outd[i]+=ind[i]*snrd[i]/ctfd[i];
00814 outd[i+1]+=ind[i+1]*snrd[i]/ctfd[i];
00815 }
00816
00817 if (nimg==1) {
00818 snrsum=snr->copy_head();
00819 float *ssnrd=snrsum->get_data();
00820
00821 for (size_t i = 0; i < sz; i+=2) { ssnrd[i]=0.0; ssnrd[i+1]=0.0; }
00822 }
00823 snr->process_inplace("math.absvalue");
00824 snrsum->add(*snr);
00825
00826 delete ctf;
00827 delete fft;
00828 delete snr;
00829 delete ctfi;
00830 }
00831
00832 EMData * CtfCAutoAverager::finish()
00833 {
00834
00835
00836
00837
00838
00839
00840
00841
00842
00843
00844 int nx=result->get_xsize();
00845 int ny=result->get_ysize();
00846 float *snrsd=snrsum->get_data();
00847 float *outd=result->get_data();
00848
00849 int rm=(ny-2)*(ny-2)/4;
00850 for (int j=0; j<ny; j++) {
00851 for (int i=0; i<nx; i+=2) {
00852 size_t ii=i+j*nx;
00853 if ((j<ny/2 && i*i/4+j*j>rm) ||(j>=ny/2 && i*i/4+(ny-j)*(ny-j)>rm) || snrsd[ii]==0) { outd[ii]=outd[ii+1]=0; continue; }
00854
00855 if (snrsd[ii]<.05) {
00856 outd[ii]*=20.0;
00857 outd[ii+1]*=20.0;
00858 }
00859 else {
00860 outd[ii]/=snrsd[ii];
00861 outd[ii+1]/=snrsd[ii];
00862 }
00863 }
00864 }
00865 result->update();
00866 result->set_attr("ptcl_repr",nimg);
00867 result->set_attr("ctf_snr_total",snrsum->calc_radial_dist(snrsum->get_ysize()/2,0,1,false));
00868 result->set_attr("ctf_wiener_filtered",0);
00869
00870
00871
00872
00873 delete snrsum;
00874 EMData *ret=result->do_ift();
00875 delete result;
00876 result=NULL;
00877 return ret;
00878 }
00879
00880 #if 0
00881 EMData *IterationAverager::average(const vector < EMData * >&image_list) const
00882 {
00883 if (image_list.size() == 0) {
00884 return 0;
00885 }
00886
00887 EMData *image0 = image_list[0];
00888
00889 int nx = image0->get_xsize();
00890 int ny = image0->get_ysize();
00891 int nz = image0->get_zsize();
00892 size_t image_size = (size_t)nx * ny * nz;
00893
00894 EMData *result = new EMData();
00895 result->set_size(nx, ny, nz);
00896
00897 EMData *sigma_image = new EMData();
00898 sigma_image->set_size(nx, ny, nz);
00899
00900 float *image0_data = image0->get_data();
00901 float *result_data = result->get_data();
00902 float *sigma_image_data = sigma_image->get_data();
00903
00904 memcpy(result_data, image0_data, image_size * sizeof(float));
00905 memcpy(sigma_image_data, image0_data, image_size * sizeof(float));
00906
00907 for (size_t j = 0; j < image_size; ++j) {
00908 sigma_image_data[j] *= sigma_image_data[j];
00909 }
00910
00911 image0->update();
00912
00913 int nc = 1;
00914 for (size_t i = 1; i < image_list.size(); i++) {
00915 EMData *image = image_list[i];
00916
00917 if (EMUtil::is_same_size(image, result)) {
00918 float *image_data = image->get_data();
00919
00920 for (size_t j = 0; j < image_size; j++) {
00921 result_data[j] += image_data[j];
00922 sigma_image_data[j] += image_data[j] * image_data[j];
00923 }
00924
00925 image->update();
00926 nc++;
00927 }
00928 }
00929
00930 float c = static_cast < float >(nc);
00931 for (size_t j = 0; j < image_size; ++j) {
00932 float f1 = sigma_image_data[j] / c;
00933 float f2 = result_data[j] / c;
00934 sigma_image_data[j] = sqrt(f1 - f2 * f2) / sqrt(c);
00935 }
00936
00937
00938 for (size_t j = 0; j < image_size; ++j) {
00939 result_data[j] /= c;
00940 }
00941
00942 result->update();
00943 sigma_image->update();
00944
00945 result->append_image("iter.hed");
00946
00947 float sigma = sigma_image->get_attr("sigma");
00948 float *sigma_image_data2 = sigma_image->get_data();
00949 float *result_data2 = result->get_data();
00950 float *d2 = new float[nx * ny];
00951 size_t sec_size = nx * ny * sizeof(float);
00952
00953 memcpy(d2, result_data2, sec_size);
00954 memcpy(sigma_image_data2, result_data2, sec_size);
00955
00956 printf("Iter sigma=%f\n", sigma);
00957
00958 for (int k = 0; k < 1000; k++) {
00959 for (int i = 1; i < nx - 1; i++) {
00960 for (int j = 1; j < ny - 1; j++) {
00961 int l = i + j * nx;
00962 float c1 = (d2[l - 1] + d2[l + 1] + d2[l - nx] + d2[l + nx]) / 4.0f - d2[l];
00963 float c2 = fabs(result_data2[l] - sigma_image_data2[l]) / sigma;
00964 result_data2[l] += c1 * Util::eman_erfc(c2) / 100.0f;
00965 }
00966 }
00967
00968 memcpy(d2, result_data2, sec_size);
00969 }
00970
00971 if( d2 )
00972 {
00973 delete[]d2;
00974 d2 = 0;
00975 }
00976
00977 sigma_image->update();
00978 if( sigma_image )
00979 {
00980 delete sigma_image;
00981 sigma_image = 0;
00982 }
00983
00984 result->update();
00985 result->append_image("iter.hed");
00986
00987 return result;
00988 }
00989 #endif
00990
00991
00992 CtfAverager::CtfAverager() :
00993 sf(0), curves(0), need_snr(false), outfile(0),
00994 image0_fft(0), image0_copy(0), snri(0), snrn(0),
00995 tdr(0), tdi(0), tn(0),
00996 filter(0), nimg(0), nx(0), ny(0), nz(0)
00997 {
00998
00999 }
01000
01001 void CtfAverager::add_image(EMData * image)
01002 {
01003 if (!image) {
01004 return;
01005 }
01006
01007 if (nimg >= 1 && !EMUtil::is_same_size(image, result)) {
01008 LOGERR("%sAverager can only process same-size Image",
01009 get_name().c_str());
01010 return;
01011 }
01012
01013 if (image->get_zsize() != 1) {
01014 LOGERR("%sAverager: Only 2D images are currently supported",
01015 get_name().c_str());
01016 }
01017
01018 string alg_name = get_name();
01019
01020 if (alg_name == "CtfCW" || alg_name == "CtfCWauto") {
01021 if (image->get_ctf() != 0 && !image->has_ctff()) {
01022 LOGERR("%sAverager: Attempted CTF Correction with no ctf parameters",
01023 get_name().c_str());
01024 }
01025 }
01026 else {
01027 if (image->get_ctf() != 0) {
01028 LOGERR("%sAverager: Attempted CTF Correction with no ctf parameters",
01029 get_name().c_str());
01030 }
01031 }
01032
01033 nimg++;
01034
01035
01036 if (nimg == 1) {
01037 image0_fft = image->do_fft();
01038
01039 nx = image0_fft->get_xsize();
01040 ny = image0_fft->get_ysize();
01041 nz = image0_fft->get_zsize();
01042
01043 result = image->copy_head();
01044 result->set_size(nx - 2, ny, nz);
01045
01046
01047 if (alg_name == "Weighting" && curves) {
01048 if (!sf) {
01049 LOGWARN("CTF curve in file will contain relative, not absolute SNR!");
01050 }
01051 curves->set_size(Ctf::CTFOS * ny / 2, 3, 1);
01052 curves->to_zero();
01053 }
01054
01055
01056 if (alg_name == "CtfC") {
01057 filter = params["filter"];
01058 if (filter == 0) {
01059 filter = 22.0f;
01060 }
01061 float apix_y = image->get_attr_dict().get("apix_y");
01062 float ds = 1.0f / (apix_y * ny * Ctf::CTFOS);
01063 filter = 1.0f / (filter * ds);
01064 }
01065
01066 if (alg_name == "CtfCWauto") {
01067 int nxy2 = nx * ny/2;
01068
01069 snri = new float[ny / 2];
01070 snrn = new float[ny / 2];
01071 tdr = new float[nxy2];
01072 tdi = new float[nxy2];
01073 tn = new float[nxy2];
01074
01075 for (int i = 0; i < ny / 2; i++) {
01076 snri[i] = 0;
01077 snrn[i] = 0;
01078 }
01079
01080 for (int i = 0; i < nxy2; i++) {
01081 tdr[i] = 1;
01082 tdi[i] = 1;
01083 tn[i] = 1;
01084 }
01085 }
01086
01087 image0_copy = image0_fft->copy_head();
01088 image0_copy->ap2ri();
01089 image0_copy->to_zero();
01090 }
01091
01092 Ctf::CtfType curve_type = Ctf::CTF_AMP;
01093 if (alg_name == "CtfCWauto") {
01094 curve_type = Ctf::CTF_AMP;
01095 }
01096
01097 float *src = image->get_data();
01098 image->ap2ri();
01099 Ctf *image_ctf = image->get_ctf();
01100 int ny2 = image->get_ysize();
01101
01102 vector<float> ctf1 = image_ctf->compute_1d(ny2,1.0f/(image_ctf->apix*image->get_ysize()), curve_type);
01103
01104 if (ctf1.size() == 0) {
01105 LOGERR("Unexpected CTF correction problem");
01106 }
01107
01108 ctf.push_back(ctf1);
01109
01110 vector<float> ctfn1;
01111 if (sf) {
01112 ctfn1 = image_ctf->compute_1d(ny2,1.0f/(image_ctf->apix*image->get_ysize()), Ctf::CTF_SNR, sf);
01113 }
01114 else {
01115 ctfn1 = image_ctf->compute_1d(ny2,1.0f/(image_ctf->apix*image->get_ysize()), Ctf::CTF_SNR);
01116 }
01117
01118 ctfn.push_back(ctfn1);
01119
01120 if (alg_name == "CtfCWauto") {
01121 int j = 0;
01122 for (int y = 0; y < ny; y++) {
01123 for (int x = 0; x < nx / 2; x++, j += 2) {
01124 #ifdef _WIN32
01125 float r = (float)_hypot((float)x, (float)(y - ny / 2.0f));
01126 #else
01127 float r = (float)hypot((float)x, (float)(y - ny / 2.0f));
01128 #endif //_WIN32
01129 int l = static_cast < int >(Util::fast_floor(r));
01130
01131 if (l >= 0 && l < ny / 2) {
01132 int k = y*nx/2 + x;
01133 tdr[k] *= src[j];
01134 tdi[k] *= src[j + 1];
01135 #ifdef _WIN32
01136 tn[k] *= (float)_hypot(src[j], src[j + 1]);
01137 #else
01138 tn[k] *= (float)hypot(src[j], src[j + 1]);
01139 #endif //_WIN32
01140 }
01141 }
01142 }
01143 }
01144
01145
01146 float *tmp_data = image0_copy->get_data();
01147
01148 int j = 0;
01149 for (int y = 0; y < ny; y++) {
01150 for (int x = 0; x < nx / 2; x++, j += 2) {
01151 float r = Ctf::CTFOS * sqrt(x * x + (y - ny / 2.0f) * (y - ny / 2.0f));
01152 int l = static_cast < int >(Util::fast_floor(r));
01153 r -= l;
01154
01155 float f = 0;
01156 if (l <= Ctf::CTFOS * ny / 2 - 2) {
01157 f = (ctf1[l] * (1 - r) + ctf1[l + 1] * r);
01158 }
01159 tmp_data[j] += src[j] * f;
01160 tmp_data[j + 1] += src[j + 1] * f;
01161 }
01162 }
01163
01164 EMData *image_fft = image->do_fft();
01165 image_fft->update();
01166 if(image_ctf) {delete image_ctf; image_ctf=0;}
01167 }
01168
01169 EMData * CtfAverager::finish()
01170 {
01171 int j = 0;
01172 for (int y = 0; y < ny; y++) {
01173 for (int x = 0; x < nx / 2; x++, j += 2) {
01174 #ifdef _WIN32
01175 float r = (float) _hypot(x, y - ny / 2.0f);
01176 #else
01177 float r = (float) hypot(x, y - ny / 2.0f);
01178 #endif
01179 int l = static_cast < int >(Util::fast_floor(r));
01180 if (l >= 0 && l < ny / 2) {
01181 int k = y*nx/2 + x;
01182 snri[l] += (tdr[k] + tdi[k]/tn[k]);
01183 snrn[l] += 1;
01184 }
01185 }
01186 }
01187
01188 for (int i = 0; i < ny / 2; i++) {
01189 snri[i] *= nimg / snrn[i];
01190 }
01191
01192 if(strcmp(outfile, "") != 0) {
01193 Util::save_data(0, 1, snri, ny / 2, outfile);
01194 }
01195
01196
01197 float *cd = 0;
01198 if (curves) {
01199 cd = curves->get_data();
01200 }
01201
01202 for (int i = 0; i < Ctf::CTFOS * ny / 2; i++) {
01203 float ctf0 = 0;
01204 for (int j = 0; j < nimg; j++) {
01205 ctf0 += ctfn[j][i];
01206 if (ctf[j][i] == 0) {
01207 ctf[j][i] = 1.0e-12f;
01208 }
01209
01210 if (curves) {
01211 cd[i] += ctf[j][i] * ctfn[j][i];
01212 cd[i + Ctf::CTFOS * ny / 2] += ctfn[j][i];
01213 cd[i + 2 * Ctf::CTFOS * ny / 2] += ctfn[j][i];
01214 }
01215 }
01216
01217 string alg_name = get_name();
01218
01219 if (alg_name == "CtfCW" && need_snr) {
01220 snr[i] = ctf0;
01221 }
01222
01223 float ctf1 = ctf0;
01224 if (alg_name == "CtfCWauto") {
01225 ctf1 = snri[i / Ctf::CTFOS];
01226 }
01227
01228 if (ctf1 <= 0.0001f) {
01229 ctf1 = 0.1f;
01230 }
01231
01232 if (alg_name == "CtfC") {
01233 for (int j = 0; j < nimg; j++) {
01234 ctf[j][i] = exp(-i * i / (filter * filter)) * ctfn[j][i] / (fabs(ctf[j][i]) * ctf1);
01235 }
01236 }
01237 else if (alg_name == "Weighting") {
01238 for (int j = 0; j < nimg; j++) {
01239 ctf[j][i] = ctfn[j][i] / ctf1;
01240 }
01241 }
01242 else if (alg_name == "CtfCW") {
01243 for (int j = 0; j < nimg; j++) {
01244 ctf[j][i] = (ctf1 / (ctf1 + 1)) * ctfn[j][i] / (ctf[j][i] * ctf1);
01245 }
01246 }
01247 else if (alg_name == "CtfCWauto") {
01248 for (int j = 0; j < nimg; j++) {
01249 ctf[j][i] = ctf1 * ctfn[j][i] / (fabs(ctf[j][i]) * ctf0);
01250 }
01251 }
01252 }
01253
01254
01255 if (curves) {
01256 for (int i = 0; i < Ctf::CTFOS * ny / 2; i++) {
01257 cd[i] /= cd[i + Ctf::CTFOS * ny / 2];
01258 }
01259 curves->update();
01260 }
01261
01262 image0_copy->update();
01263
01264 float *result_data = result->get_data();
01265 EMData *tmp_ift = image0_copy->do_ift();
01266 float *tmp_ift_data = tmp_ift->get_data();
01267 memcpy(result_data, tmp_ift_data, (nx - 2) * ny * sizeof(float));
01268
01269 tmp_ift->update();
01270 result->update();
01271 result->set_attr("ptcl_repr",nimg);
01272
01273 if( image0_copy )
01274 {
01275 delete image0_copy;
01276 image0_copy = 0;
01277 }
01278
01279 if (snri) {
01280 delete[]snri;
01281 snri = 0;
01282 }
01283
01284 if (snrn) {
01285 delete[]snrn;
01286 snrn = 0;
01287 }
01288
01289 if( snri )
01290 {
01291 delete [] snri;
01292 snri = 0;
01293 }
01294 if( snrn )
01295 {
01296 delete [] snrn;
01297 snrn = 0;
01298 }
01299 if( tdr )
01300 {
01301 delete [] tdr;
01302 tdr = 0;
01303 }
01304 if( tdi )
01305 {
01306 delete [] tdi;
01307 tdi = 0;
01308 }
01309 if( tn )
01310 {
01311 delete [] tn;
01312 tn = 0;
01313 }
01314
01315 return result;
01316 }
01317
01318 #if 0
01319 EMData *CtfAverager::average(const vector < EMData * >&image_list) const
01320 {
01321 if (image_list.size() == 0) {
01322 return 0;
01323 }
01324
01325 EMData *image0 = image_list[0];
01326 if (image0->get_zsize() != 1) {
01327 LOGERR("Only 2D images are currently supported");
01328 return 0;
01329 }
01330
01331 string alg_name = get_name();
01332
01333 if (alg_name == "CtfCW" || alg_name == "CtfCWauto") {
01334 if (image0->get_ctf() != 0 && !image0->has_ctff()) {
01335 LOGERR("Attempted CTF Correction with no ctf parameters");
01336 return 0;
01337 }
01338 }
01339 else {
01340 if (image0->get_ctf() != 0) {
01341 LOGERR("Attempted CTF Correction with no ctf parameters");
01342 return 0;
01343 }
01344 }
01345
01346 size_t num_images = image_list.size();
01347 vector < float >*ctf = new vector < float >[num_images];
01348 vector < float >*ctfn = new vector < float >[num_images];
01349 float **src = new float *[num_images];
01350
01351 Ctf::CtfType curve_type = Ctf::CTF_ABS_AMP;
01352 if (alg_name == "CtfCWauto") {
01353 curve_type = Ctf::CTF_AMP;
01354 }
01355
01356 for (size_t i = 0; i < num_images; i++) {
01357 EMData *image = image_list[i]->do_fft();
01358 image->ap2ri();
01359 src[i] = image->get_data();
01360 Ctf *image_ctf = image->get_ctf();
01361 int ny = image->get_ysize();
01362 ctf[i] = image_ctf->compute_1d(ny, curve_type);
01363
01364 if (ctf[i].size() == 0) {
01365 LOGERR("Unexpected CTF correction problem");
01366 return 0;
01367 }
01368
01369 if (sf) {
01370 ctfn[i] = image_ctf->compute_1d(ny, Ctf::CTF_ABS_SNR, sf);
01371 }
01372 else { result->set_attr("ptcl_repr",nimg);
01373 ctfn[i] = image_ctf->compute_1d(ny, Ctf::CTF_RELATIVE_SNR);
01374 }
01375
01376 if(image_ctf) {delete image_ctf; image_ctf=0;}
01377 }
01378
01379 EMData *image0_fft = image0->do_fft();
01380
01381 int nx = image0_fft->get_xsize();
01382 int ny = image0_fft->get_ysize();
01383 int nz = image0_fft->get_zsize();
01384
01385 EMData *result = new EMData();
01386 result->set_size(nx - 2, ny, nz);
01387
01388 float *cd = 0;
01389 if (alg_name == "Weighting" && curves) {
01390 if (!sf) {
01391 LOGWARN("CTF curve in file will contain relative, not absolute SNR!");
01392 }
01393 curves->set_size(Ctf::CTFOS * ny / 2, 3, 1);
01394 curves->to_zero();
01395 cd = curves->get_data();
01396 }
01397
01398 float filter = 0;
01399 if (alg_name == "CtfC") {
01400 filter = params["filter"];
01401 if (filter == 0) {
01402 filter = 22.0f;
01403 }
01404 float apix_y = image0->get_attr_dict().get("apix_y");
01405 float ds = 1.0f / (apix_y * ny * Ctf::CTFOS);
01406 filter = 1.0f / (filter * ds);
01407 }
01408
01409 float *snri = 0;
01410 float *snrn = 0;
01411
01412 if (alg_name == "CtfCWauto") {
01413 snri = new float[ny / 2];
01414 snrn = new float[ny / 2];
01415
01416 for (int i = 0; i < ny / 2; i++) {
01417 snri[i] = 0;
01418 snrn[i] = 0;
01419 }
01420
01421 int j = 0;
01422 for (int y = 0; y < ny; y++) {
01423 for (int x = 0; x < nx / 2; x++, j += 2) {
01424 float r = hypot(x, y - ny / 2.0f);
01425 int l = static_cast < int >(Util::fast_floor(r));
01426
01427 if (l >= 0 && l < ny / 2) {
01428 float tdr = 1;
01429 float tdi = 1;
01430 float tn = 1;
01431
01432 for (size_t i = 0; i < num_images; i++) {
01433 tdr *= src[i][j];
01434 tdi *= src[i][j + 1];
01435 tn *= hypot(src[i][j], src[i][j + 1]);
01436 }
01437
01438 tdr += tdi / tn;
01439 snri[l] += tdr;
01440 snrn[l] += 1;
01441 }
01442 }
01443 }
01444
01445 for (int i = 0; i < ny / 2; i++) {
01446 snri[i] *= num_images / snrn[i];
01447 }
01448 if (outfile != "") {
01449 Util::save_data(0, 1, snri, ny / 2, outfile);
01450 }
01451 }
01452
01453 for (int i = 0; i < Ctf::CTFOS * ny / 2; i++) {
01454 float ctf0 = 0;
01455 for (size_t j = 0; j < num_images; j++) {
01456 ctf0 += ctfn[j][i];
01457 if (ctf[j][i] == 0) {
01458 ctf[j][i] = 1.0e-12;
01459 }
01460
01461 if (curves) {
01462 cd[i] += ctf[j][i] * ctfn[j][i];
01463 cd[i + Ctf::CTFOS * ny / 2] += ctfn[j][i];
01464 cd[i + 2 * Ctf::CTFOS * ny / 2] += ctfn[j][i];
01465 }
01466 }
01467
01468 if (alg_name == "CtfCW" && need_snr) {
01469 snr[i] = ctf0;
01470 }
01471
01472 float ctf1 = ctf0;
01473 if (alg_name == "CtfCWauto") {
01474 ctf1 = snri[i / Ctf::CTFOS];
01475 }
01476
01477 if (ctf1 <= 0.0001f) {
01478 ctf1 = 0.1f;
01479 }
01480
01481 if (alg_name == "CtfC") {
01482 for (size_t j = 0; j < num_images; j++) {
01483 ctf[j][i] = exp(-i * i / (filter * filter)) * ctfn[j][i] / (fabs(ctf[j][i]) * ctf1);
01484 }
01485 }
01486 else if (alg_name == "Weighting") {
01487 for (size_t j = 0; j < num_images; j++) {
01488 ctf[j][i] = ctfn[j][i] / ctf1;
01489 }
01490 }
01491 else if (alg_name == "CtfCW") {
01492 for (size_t j = 0; j < num_images; j++) {
01493 ctf[j][i] = (ctf1 / (ctf1 + 1)) * ctfn[j][i] / (ctf[j][i] * ctf1);
01494 }
01495 }
01496 else if (alg_name == "CtfCWauto") {
01497 for (size_t j = 0; j < num_images; j++) {
01498 ctf[j][i] = ctf1 * ctfn[j][i] / (fabs(ctf[j][i]) * ctf0);
01499 }
01500 }
01501 }
01502
01503
01504 if (curves) {
01505 for (int i = 0; i < Ctf::CTFOS * ny / 2; i++) {
01506 cd[i] /= cd[i + Ctf::CTFOS * ny / 2];
01507 }
01508 curves->update();
01509 }
01510
01511 EMData *image0_copy = image0_fft->copy_head();
01512 image0_copy->ap2ri();
01513
01514 float *tmp_data = image0_copy->get_data();
01515
01516 int j = 0;
01517 for (int y = 0; y < ny; y++) {
01518 for (int x = 0; x < nx / 2; x++, j += 2) {
01519 float r = Ctf::CTFOS * sqrt(x * x + (y - ny / 2.0f) * (y - ny / 2.0f));
01520 int l = static_cast < int >(Util::fast_floor(r));
01521 r -= l;
01522
01523 tmp_data[j] = 0;
01524 tmp_data[j + 1] = 0;
01525
01526 for (size_t i = 0; i < num_images; i++) {
01527 float f = 0;
01528 if (l <= Ctf::CTFOS * ny / 2 - 2) {
01529 f = (ctf[i][l] * (1 - r) + ctf[i][l + 1] * r);
01530 }
01531 tmp_data[j] += src[i][j] * f;
01532 tmp_data[j + 1] += src[i][j + 1] * f;
01533 }
01534 }
01535 }
01536
01537 image0_copy->update();
01538
01539 float *result_data = result->get_data();
01540 EMData *tmp_ift = image0_copy->do_ift();
01541 float *tmp_ift_data = tmp_ift->get_data();
01542 memcpy(result_data, tmp_ift_data, (nx - 2) * ny * sizeof(float));
01543
01544 tmp_ift->update();
01545
01546 if( image0_copy )
01547 {
01548 delete image0_copy;
01549 image0_copy = 0;
01550 }
01551
01552 for (size_t i = 0; i < num_images; i++) {
01553 EMData *img = image_list[i]->do_fft();
01554 img->update();
01555 }
01556
01557 if( src )
01558 {
01559 delete[]src;
01560 src = 0;
01561 }
01562
01563 if( ctf )
01564 {
01565 delete[]ctf;
01566 ctf = 0;
01567 }
01568
01569 if( ctfn )
01570 {
01571 delete[]ctfn;
01572 ctfn = 0;
01573 }
01574
01575 if (snri) {
01576 delete[]snri;
01577 snri = 0;
01578 }
01579
01580 if (snrn) {
01581 delete[]snrn;
01582 snrn = 0;
01583 }
01584
01585 result->update();
01586 return result;
01587 }
01588 #endif
01589
01590
01591 void EMAN::dump_averagers()
01592 {
01593 dump_factory < Averager > ();
01594 }
01595
01596 map<string, vector<string> > EMAN::dump_averagers_list()
01597 {
01598 return dump_factory_list < Averager > ();
01599 }