#include <util.h>
Inheritance diagram for EMAN::Util::KaiserBessel:
Public Member Functions | |
KaiserBessel (float alpha_, int K, float r_, float v_, int N_, float vtable_=0.f, int ntable_=5999) | |
virtual | ~KaiserBessel () |
float | I0table_maxerror () |
Compute the maximum error in the table. | |
vector< float > | dump_table () |
virtual float | sinhwin (float x) const |
Kaiser-Bessel Sinh window function. | |
virtual float | i0win (float x) const |
Kaiser-Bessel I0 window function. | |
float | i0win_tab (float x) const |
Kaiser-Bessel I0 window function (uses table lookup). | |
int | get_window_size () const |
Return the size of the I0 window. | |
kbsinh_win | get_kbsinh_win () |
Sinh window function object factory. | |
kbi0_win | get_kbi0_win () |
I0 window function object factory. | |
Protected Member Functions | |
virtual void | build_I0table () |
2*pi*alpha*r*vadjust | |
Protected Attributes | |
float | alpha |
float | v |
float | r |
int | N |
Kaiser-Bessel parameters. | |
int | K |
size in Ix-space | |
float | vtable |
I0 window size. | |
int | ntable |
table I0 non-zero domain maximum | |
vector< float > | i0table |
float | dtable |
float | alphar |
table spacing | |
float | fac |
alpha*r | |
float | vadjust |
2*pi*alpha*r*v | |
float | facadj |
float | fltb |
Tabulate I0 window for speed. | |
Classes | |
class | kbi0_win |
I0 window function object. More... | |
class | kbsinh_win |
Sinh window function object. More... |
(It's a class so that the windowing parameters may be instantiated and held in the instance object.)
The I0 version can be tabulated and interpolated upon demand, but the max error needs to be checked. The "vtable" parameter corresponds to the maximum value of x for which the I0 window is non-zero. Setting "vtable" different from "v" corresponds to a change in units of x. In practice, it is often handy to replace x in some sort of absolute units with x described in terms of grid intervals.
The get_kbsinh_win and get_kbi0_win functions return single-argument function objects, which is what a generic routine is likely to want.
Definition at line 307 of file util.h.
Util::KaiserBessel::KaiserBessel | ( | float | alpha_, | |
int | K, | |||
float | r_, | |||
float | v_, | |||
int | N_, | |||
float | vtable_ = 0.f , |
|||
int | ntable_ = 5999 | |||
) |
Definition at line 1982 of file util_sparx.cpp.
References alpha, alphar, build_I0table(), fac, facadj, K, r, twopi, v, vadjust, and vtable.
01984 : alpha(alpha_), v(v_), r(r_), N(N_), K(K_), vtable(vtable_), 01985 ntable(ntable_) { 01986 // Default values are alpha=1.25, K=6, r=0.5, v = K/2 01987 if (0.f == v) v = float(K)/2; 01988 if (0.f == vtable) vtable = v; 01989 alphar = alpha*r; 01990 fac = static_cast<float>(twopi)*alphar*v; 01991 vadjust = 1.0f*v; 01992 facadj = static_cast<float>(twopi)*alphar*vadjust; 01993 build_I0table(); 01994 }
virtual EMAN::Util::KaiserBessel::~KaiserBessel | ( | ) | [inline, virtual] |
void Util::KaiserBessel::build_I0table | ( | ) | [protected, virtual] |
2*pi*alpha*r*vadjust
Reimplemented in EMAN::Util::FakeKaiserBessel.
Definition at line 2005 of file util_sparx.cpp.
References facadj, fltb, i0table, K, N, ntable, EMAN::Util::round(), sqrt(), and vadjust.
Referenced by KaiserBessel().
02005 { 02006 i0table.resize(ntable+1); // i0table[0:ntable] 02007 int ltab = int(round(float(ntable)/1.25f)); 02008 fltb = float(ltab)/(K/2); 02009 float val0 = static_cast<float>(gsl_sf_bessel_I0(facadj)); 02010 for (int i=ltab+1; i <= ntable; i++) i0table[i] = 0.f; 02011 for (int i=0; i <= ltab; i++) { 02012 float s = float(i)/fltb/N; 02013 if (s < vadjust) { 02014 float rt = sqrt(1.f - pow(s/vadjust, 2)); 02015 i0table[i] = static_cast<float>(gsl_sf_bessel_I0(facadj*rt))/val0; 02016 } else { 02017 i0table[i] = 0.f; 02018 } 02019 // cout << " "<<s*N<<" "<<i0table[i] <<endl; 02020 } 02021 }
vector<float> EMAN::Util::KaiserBessel::dump_table | ( | ) | [inline] |
kbi0_win EMAN::Util::KaiserBessel::get_kbi0_win | ( | ) | [inline] |
kbsinh_win EMAN::Util::KaiserBessel::get_kbsinh_win | ( | ) | [inline] |
int EMAN::Util::KaiserBessel::get_window_size | ( | ) | const [inline] |
Return the size of the I0 window.
Definition at line 351 of file util.h.
Referenced by EMAN::EMData::extract_plane(), EMAN::EMData::extractline(), EMAN::EMData::extractpoint(), EMAN::EMData::get_pixel_conv(), EMAN::Util::get_pixel_conv_new(), EMAN::Util::get_pixel_conv_new_background(), EMAN::EMData::rot_scale_conv(), EMAN::EMData::rot_scale_conv7(), EMAN::EMData::rot_scale_conv_new(), and EMAN::EMData::rot_scale_conv_new_background().
float Util::KaiserBessel::I0table_maxerror | ( | ) |
Compute the maximum error in the table.
Definition at line 2023 of file util_sparx.cpp.
References i0table, and ntable.
02023 { 02024 float maxdiff = 0.f; 02025 for (int i = 1; i <= ntable; i++) { 02026 float diff = fabs(i0table[i] - i0table[i-1]); 02027 if (diff > maxdiff) maxdiff = diff; 02028 } 02029 return maxdiff; 02030 }
float Util::KaiserBessel::i0win | ( | float | x | ) | const [virtual] |
Kaiser-Bessel I0 window function.
Reimplemented in EMAN::Util::FakeKaiserBessel.
Definition at line 1996 of file util_sparx.cpp.
References facadj, sqrt(), and vadjust.
Referenced by EMAN::Processor::EMFourierFilterFunc().
01996 { 01997 float val0 = float(gsl_sf_bessel_I0(facadj)); 01998 float absx = fabs(x); 01999 if (absx > vadjust) return 0.f; 02000 float rt = sqrt(1.f - pow(absx/vadjust, 2)); 02001 float res = static_cast<float>(gsl_sf_bessel_I0(facadj*rt))/val0; 02002 return res; 02003 }
float EMAN::Util::KaiserBessel::i0win_tab | ( | float | x | ) | const [inline] |
Kaiser-Bessel I0 window function (uses table lookup).
Definition at line 338 of file util.h.
Referenced by EMAN::EMData::extract_plane(), EMAN::EMData::extractline(), EMAN::EMData::extractpoint(), EMAN::Util::extractpoint2(), EMAN::EMData::get_pixel_conv(), EMAN::Util::get_pixel_conv_new(), EMAN::Util::get_pixel_conv_new_background(), EMAN::EMData::rot_scale_conv(), and EMAN::EMData::rot_scale_conv7().
float Util::KaiserBessel::sinhwin | ( | float | x | ) | const [virtual] |
Kaiser-Bessel Sinh window function.
Reimplemented in EMAN::Util::FakeKaiserBessel.
Definition at line 2032 of file util_sparx.cpp.
References alphar, fac, and sqrt().
Referenced by EMAN::EMData::divkbsinh(), and EMAN::Processor::EMFourierFilterFunc().
02032 { 02033 float val0 = sinh(fac)/fac; 02034 float absx = fabs(x); 02035 if (0.0 == x) { 02036 float res = 1.0f; 02037 return res; 02038 } else if (absx == alphar) { 02039 return 1.0f/val0; 02040 } else if (absx < alphar) { 02041 float rt = sqrt(1.0f - pow((x/alphar), 2)); 02042 float facrt = fac*rt; 02043 float res = (sinh(facrt)/facrt)/val0; 02044 return res; 02045 } else { 02046 float rt = sqrt(pow((x/alphar),2) - 1.f); 02047 float facrt = fac*rt; 02048 float res = (sin(facrt)/facrt)/val0; 02049 return res; 02050 } 02051 }
float EMAN::Util::KaiserBessel::alpha [protected] |
float EMAN::Util::KaiserBessel::alphar [protected] |
table spacing
Definition at line 317 of file util.h.
Referenced by KaiserBessel(), EMAN::Util::FakeKaiserBessel::sinhwin(), and sinhwin().
float EMAN::Util::KaiserBessel::dtable [protected] |
float EMAN::Util::KaiserBessel::fac [protected] |
alpha*r
Definition at line 318 of file util.h.
Referenced by KaiserBessel(), EMAN::Util::FakeKaiserBessel::sinhwin(), and sinhwin().
float EMAN::Util::KaiserBessel::facadj [protected] |
Definition at line 320 of file util.h.
Referenced by EMAN::Util::FakeKaiserBessel::build_I0table(), build_I0table(), EMAN::Util::FakeKaiserBessel::i0win(), i0win(), and KaiserBessel().
float EMAN::Util::KaiserBessel::fltb [protected] |
Tabulate I0 window for speed.
Definition at line 322 of file util.h.
Referenced by EMAN::Util::FakeKaiserBessel::build_I0table(), and build_I0table().
vector<float> EMAN::Util::KaiserBessel::i0table [protected] |
Definition at line 315 of file util.h.
Referenced by EMAN::Util::FakeKaiserBessel::build_I0table(), build_I0table(), and I0table_maxerror().
int EMAN::Util::KaiserBessel::K [protected] |
size in Ix-space
Definition at line 312 of file util.h.
Referenced by EMAN::Util::FakeKaiserBessel::build_I0table(), build_I0table(), and KaiserBessel().
int EMAN::Util::KaiserBessel::N [protected] |
Kaiser-Bessel parameters.
Definition at line 311 of file util.h.
Referenced by EMAN::Util::FakeKaiserBessel::build_I0table(), and build_I0table().
int EMAN::Util::KaiserBessel::ntable [protected] |
table I0 non-zero domain maximum
Definition at line 314 of file util.h.
Referenced by EMAN::Util::FakeKaiserBessel::build_I0table(), build_I0table(), and I0table_maxerror().
float EMAN::Util::KaiserBessel::r [protected] |
float EMAN::Util::KaiserBessel::v [protected] |
float EMAN::Util::KaiserBessel::vadjust [protected] |
2*pi*alpha*r*v
Definition at line 319 of file util.h.
Referenced by EMAN::Util::FakeKaiserBessel::build_I0table(), build_I0table(), EMAN::Util::FakeKaiserBessel::i0win(), i0win(), and KaiserBessel().
float EMAN::Util::KaiserBessel::vtable [protected] |