#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. |
(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 305 of file util.h.
|
Definition at line 1986 of file util_sparx.cpp. References alpha, alphar, build_I0table(), fac, facadj, K, v, v, vadjust, and vtable. 01988 : alpha(alpha_), v(v_), r(r_), N(N_), K(K_), vtable(vtable_), 01989 ntable(ntable_) { 01990 // Default values are alpha=1.25, K=6, r=0.5, v = K/2 01991 if (0.f == v) v = float(K)/2; 01992 if (0.f == vtable) vtable = v; 01993 alphar = alpha*r; 01994 fac = static_cast<float>(twopi)*alphar*v; 01995 vadjust = 1.0f*v; 01996 facadj = static_cast<float>(twopi)*alphar*vadjust; 01997 build_I0table(); 01998 }
|
|
Definition at line 325 of file util.h. 00439 {
|
|
2*pi*alpha*r*vadjust
Reimplemented in EMAN::Util::FakeKaiserBessel. Definition at line 2009 of file util_sparx.cpp. References facadj, fltb, i0table, K, ntable, EMAN::Util::round(), sqrt(), and vadjust. Referenced by KaiserBessel(). 02009 { 02010 i0table.resize(ntable+1); // i0table[0:ntable] 02011 int ltab = int(round(float(ntable)/1.25f)); 02012 fltb = float(ltab)/(K/2); 02013 float val0 = static_cast<float>(gsl_sf_bessel_I0(facadj)); 02014 for (int i=ltab+1; i <= ntable; i++) i0table[i] = 0.f; 02015 for (int i=0; i <= ltab; i++) { 02016 float s = float(i)/fltb/N; 02017 if (s < vadjust) { 02018 float rt = sqrt(1.f - pow(s/vadjust, 2)); 02019 i0table[i] = static_cast<float>(gsl_sf_bessel_I0(facadj*rt))/val0; 02020 } else { 02021 i0table[i] = 0.f; 02022 } 02023 // cout << " "<<s*N<<" "<<i0table[i] <<endl; 02024 } 02025 }
|
|
Definition at line 328 of file util.h. References EMAN::Vec3f. 00439 {
|
|
I0 window function object factory.
Definition at line 375 of file util.h. References x. 00439 {
|
|
Sinh window function object factory.
Definition at line 361 of file util.h. 00439 {
|
|
Return the size of the I0 window.
Definition at line 349 of file util.h. Referenced by EMAN::EMData::extract_plane(), EMAN::EMData::extract_plane_rect(), EMAN::EMData::extract_plane_rect_fast(), 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(), and EMAN::EMData::rot_scale_conv7(). 00439 {
|
|
Compute the maximum error in the table.
Definition at line 2027 of file util_sparx.cpp. References i0table. 02027 { 02028 float maxdiff = 0.f; 02029 for (int i = 1; i <= ntable; i++) { 02030 float diff = fabs(i0table[i] - i0table[i-1]); 02031 if (diff > maxdiff) maxdiff = diff; 02032 } 02033 return maxdiff; 02034 }
|
|
Kaiser-Bessel I0 window function.
Reimplemented in EMAN::Util::FakeKaiserBessel. Definition at line 2000 of file util_sparx.cpp. References facadj, sqrt(), vadjust, and x. Referenced by EMAN::Processor::EMFourierFilterFunc(). 02000 { 02001 float val0 = float(gsl_sf_bessel_I0(facadj)); 02002 float absx = fabs(x); 02003 if (absx > vadjust) return 0.f; 02004 float rt = sqrt(1.f - pow(absx/vadjust, 2)); 02005 float res = static_cast<float>(gsl_sf_bessel_I0(facadj*rt))/val0; 02006 return res; 02007 }
|
|
Kaiser-Bessel I0 window function (uses table lookup).
Definition at line 336 of file util.h. Referenced by EMAN::EMData::extract_plane(), EMAN::EMData::extract_plane_rect(), EMAN::EMData::extract_plane_rect_fast(), 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(). 00439 {
|
|
Kaiser-Bessel Sinh window function.
Reimplemented in EMAN::Util::FakeKaiserBessel. Definition at line 2036 of file util_sparx.cpp. References alphar, fac, sqrt(), and x. Referenced by EMAN::EMData::divkbsinh(), EMAN::EMData::divkbsinh_rect(), and EMAN::Processor::EMFourierFilterFunc(). 02036 { 02037 float val0 = sinh(fac)/fac; 02038 float absx = fabs(x); 02039 if (0.0 == x) { 02040 float res = 1.0f; 02041 return res; 02042 } else if (absx == alphar) { 02043 return 1.0f/val0; 02044 } else if (absx < alphar) { 02045 float rt = sqrt(1.0f - pow((x/alphar), 2)); 02046 float facrt = fac*rt; 02047 float res = (sinh(facrt)/facrt)/val0; 02048 return res; 02049 } else { 02050 float rt = sqrt(pow((x/alphar),2) - 1.f); 02051 float facrt = fac*rt; 02052 float res = (sin(facrt)/facrt)/val0; 02053 return res; 02054 } 02055 }
|
|
Definition at line 308 of file util.h. Referenced by KaiserBessel(). |
|
table spacing
Definition at line 315 of file util.h. Referenced by KaiserBessel(), and sinhwin(). |
|
|
|
alpha*r
Definition at line 316 of file util.h. Referenced by KaiserBessel(), and sinhwin(). |
|
Definition at line 318 of file util.h. Referenced by build_I0table(), i0win(), and KaiserBessel(). |
|
Tabulate I0 window for speed.
Definition at line 320 of file util.h. Referenced by build_I0table(). |
|
Definition at line 313 of file util.h. Referenced by build_I0table(), and I0table_maxerror(). |
|
size in Ix-space
Definition at line 310 of file util.h. Referenced by build_I0table(), and KaiserBessel(). |
|
Kaiser-Bessel parameters.
|
|
table I0 non-zero domain maximum
Definition at line 312 of file util.h. Referenced by build_I0table(). |
|
|
|
Definition at line 308 of file util.h. Referenced by KaiserBessel(). |
|
2*pi*alpha*r*v
Definition at line 317 of file util.h. Referenced by build_I0table(), i0win(), and KaiserBessel(). |
|
I0 window size.
Definition at line 311 of file util.h. Referenced by KaiserBessel(). |