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 #ifndef eman__util_h__
00037 #define eman__util_h__ 1
00038
00039 #ifdef _WIN32
00040 #pragma warning(disable:4819)
00041 #endif //_WIN32
00042
00043 #include "sparx/emconstants.h"
00044 #include "exception.h"
00045 #include <vector>
00046 #include <iostream>
00047
00048 #include <string>
00049 using std::string;
00050
00051
00052 #include <boost/multi_array.hpp>
00053 #include <boost/tuple/tuple.hpp>
00054
00055 #include "vec3.h"
00056
00057 #ifdef WIN32
00058 #include <windows.h>
00059 #define M_PI 3.14159265358979323846f
00060
00061 #endif
00062
00063 using std::string;
00064 using std::vector;
00065 using std::ostream;
00066 using std::cout;
00067 using std::endl;
00068
00069
00070 namespace EMAN
00071 {
00072 class EMData;
00073 class Dict;
00074
00075 typedef boost::multi_array<int, 3> MIArray3D;
00076
00080 class Util
00081 {
00083 #include "sparx/util_sparx.h"
00084
00085 public:
00091 static void ap2ri(float *data, size_t n);
00092
00097 static void flip_complex_phase(float *data, size_t n);
00098
00106 static void rotate_phase_origin(float *data, size_t nx, size_t ny, size_t nz);
00107
00115 static int file_lock_wait(FILE * file);
00116
00117
00123 static bool check_file_by_magic(const void *first_block, const char *magic);
00124
00128 static bool is_file_exist(const string & filename);
00129
00135 static void flip_image(float *data, size_t nx, size_t ny);
00136
00142 static vector<EMData *> svdcmp(const vector<EMData *> &data,int nvec);
00143
00144
00149 static string str_to_lower(const string& s);
00150
00159 static bool sstrncmp(const char *s1, const char *s2);
00160
00165 static string int2str(int n);
00166
00175 static string get_line_from_string(char **str);
00176
00188 static bool get_str_float(const char *str, const char *float_var, float *p_val);
00189
00190
00203 static bool get_str_float(const char *str, const char *float_var,
00204 float *p_v1, float *p_v2);
00205
00224 static bool get_str_float(const char *str, const char *float_var,
00225 int *p_nvalues, float *p_v1, float *p_v2);
00226
00238 static bool get_str_int(const char *str, const char *int_var, int *p_val);
00239
00252 static bool get_str_int(const char *str, const char *int_var,
00253 int *p_v1, int *p_v2);
00254
00273 static bool get_str_int(const char *str, const char *int_var,
00274 int *p_nvalues, int *p_v1, int *p_v2);
00275
00286 static string change_filename_ext(const string& old_filename,
00287 const string & new_ext);
00288
00295 static string remove_filename_ext(const string& filename);
00296
00302 static string get_filename_ext(const string& filename);
00303
00310 static string sbasename(const string & filename);
00311
00324 static void calc_least_square_fit(size_t nitems, const float *data_x,
00325 const float *data_y, float *p_slope,
00326 float *p_intercept, bool ignore_zero,float absmax=0);
00327
00334 static Vec3f calc_bilinear_least_square(const vector<float> &points);
00335
00346 static void save_data(const vector < float >&x_array,
00347 const vector < float >&y_array,
00348 const string & filename);
00349
00358 static void save_data(float x0, float dx,
00359 const vector < float >&y_array,
00360 const string & filename);
00361
00371 static void save_data(float x0, float dx, float *y_array,
00372 size_t array_size, const string & filename);
00373
00374
00383 static void sort_mat(float *left, float *right, int *leftPerm,
00384 int *rightPerm);
00385
00386
00390 static unsigned long int get_randnum_seed();
00391
00395 static void set_randnum_seed(unsigned long int seed);
00396
00402 static int get_irand(int low, int high);
00403
00409 static float get_frand(int low, int high);
00410
00416 static float get_frand(float low, float high);
00417
00423 static float get_frand(double low, double high);
00424
00431 static float get_gauss_rand(float mean, float sigma);
00432
00437 static inline int round(float x)
00438 {
00439 if (x < 0) {
00440 return (int) (x - 0.5f);
00441 }
00442 return (int) (x + 0.5f);
00443 }
00444
00449 static inline int round(double x)
00450 {
00451 if (x < 0) {
00452 return (int) (x - 0.5);
00453 }
00454 return (int) (x + 0.5);
00455 }
00456
00463 static inline float linear_interpolate(float p1, float p2, float t)
00464 {
00465 return (1-t) * p1 + t * p2;
00466 }
00467
00477 static inline float bilinear_interpolate(float p1, float p2, float p3,
00478 float p4, float t, float u)
00479 {
00480 return (1-t) * (1-u) * p1 + t * (1-u) * p2 + (1-t) * u * p3 + t * u * p4;
00481 }
00482
00498 static inline float trilinear_interpolate(float p1, float p2, float p3,
00499 float p4, float p5, float p6,
00500 float p7, float p8, float t,
00501 float u, float v)
00502 {
00503 return ((1 - t) * (1 - u) * (1 - v) * p1 + t * (1 - u) * (1 - v) * p2
00504 + (1 - t) * u * (1 - v) * p3 + t * u * (1 - v) * p4
00505 + (1 - t) * (1 - u) * v * p5 + t * (1 - u) * v * p6
00506 + (1 - t) * u * v * p7 + t * u * v * p8);
00507 }
00508
00515 static void find_max(const float *data, size_t nitems,
00516 float *p_max_val, int *p_max_index = 0);
00517
00528 static void find_min_and_max(const float *data, size_t nitems,
00529 float *p_max_val, float *p_min_val,
00530 int *p_max_index = 0, int *p_min_index = 0);
00531
00532
00537 static Dict get_stats( const vector<float>& data );
00538
00542 static Dict get_stats_cstyle( const vector<float>& data );
00543
00550 static int calc_best_fft_size(int low);
00551
00552
00553 static EMData* calc_bessel(const int n, const float& x);
00554
00559 static inline int square(int n)
00560 {
00561 return (n * n);
00562 }
00563
00568 static inline float square(float x)
00569 {
00570 return (x * x);
00571 }
00572
00577 static inline float square(double x)
00578 {
00579 return (float)(x * x);
00580 }
00581
00587 static inline float square_sum(float x, float y)
00588 {
00589 return (float)(x * x + y * y);
00590 }
00591
00597 static inline float hypot2(float x, float y)
00598 {
00599 return sqrtf(x * x + y * y);
00600 }
00601
00608 static inline int hypot3sq(int x, int y, int z)
00609 {
00610 return ((x * x + y * y + z * z));
00611 }
00612
00619 static inline float hypot3sq(float x, float y, float z)
00620 {
00621 return (x * x + y * y + z * z);
00622 }
00623
00630 static inline float hypot3(int x, int y, int z)
00631 {
00632 return sqrtf((float)(x * x + y * y + z * z));
00633 }
00634
00641 static inline float hypot3(float x, float y, float z)
00642 {
00643 return sqrtf(x * x + y * y + z * z);
00644 }
00645
00652 static inline float hypot3(double x, double y, double z)
00653 {
00654 return (float) sqrt(x * x + y * y + z * z);
00655 }
00656
00662 static float hypot_fast(int x, int y);
00663
00669 static short hypot_fast_int(int x, int y);
00670
00677 static inline int fast_floor(float x)
00678 {
00679 if (x < 0) {
00680 return ((int) x - 1);
00681 }
00682 return (int) x;
00683 }
00689 static float fast_exp(const float &f) ;
00690
00696 static float fast_acos(const float &f) ;
00697
00706 static inline float agauss(float a, float dx, float dy, float dz, float d)
00707 {
00708 return (a * exp(-(dx * dx + dy * dy + dz * dz) / d));
00709 }
00710
00716 static inline int get_min(int f1, int f2)
00717 {
00718 return (f1 < f2 ? f1 : f2);
00719 }
00720
00727 static inline int get_min(int f1, int f2, int f3)
00728 {
00729 if (f1 <= f2 && f1 <= f3) {
00730 return f1;
00731 }
00732 if (f2 <= f1 && f2 <= f3) {
00733 return f2;
00734 }
00735 return f3;
00736 }
00737
00743 static inline float get_min(float f1, float f2)
00744 {
00745 return (f1 < f2 ? f1 : f2);
00746 }
00747
00754 static inline float get_min(float f1, float f2, float f3)
00755 {
00756 if (f1 <= f2 && f1 <= f3) {
00757 return f1;
00758 }
00759 if (f2 <= f1 && f2 <= f3) {
00760 return f2;
00761 }
00762 return f3;
00763 }
00764
00765
00773 static inline float get_min(float f1, float f2, float f3, float f4)
00774 {
00775 float m = f1;
00776 if (f2 < m) {
00777 m = f2;
00778 }
00779 if (f3 < m) {
00780 m = f3;
00781 }
00782 if (f4 < m) {
00783 m = f4;
00784 }
00785 return m;
00786 }
00787
00793 static inline float get_max(float f1, float f2)
00794 {
00795 return (f1 < f2 ? f2 : f1);
00796 }
00797
00804 static inline float get_max(float f1, float f2, float f3)
00805 {
00806 if (f1 >= f2 && f1 >= f3) {
00807 return f1;
00808 }
00809 if (f2 >= f1 && f2 >= f3) {
00810 return f2;
00811 }
00812 return f3;
00813 }
00814
00822 static inline float get_max(float f1, float f2, float f3, float f4)
00823 {
00824 float m = f1;
00825 if (f2 > m) {
00826 m = f2;
00827 }
00828 if (f3 > m) {
00829 m = f3;
00830 }
00831 if (f4 > m) {
00832 m = f4;
00833 }
00834 return m;
00835 }
00836
00844 static inline float angle_sub_2pi(float x, float y)
00845 {
00846 float r = fmod(fabs(x - y), (float) (2 * M_PI));
00847 if (r > M_PI) {
00848 r = (float) (2.0 * M_PI - r);
00849 }
00850
00851 return r;
00852 }
00853
00861 static inline float angle_sub_pi(float x, float y)
00862 {
00863 float r = fmod(fabs(x - y), (float) M_PI);
00864 if (r > M_PI / 2.0) {
00865 r = (float)(M_PI - r);
00866 }
00867 return r;
00868 }
00869
00876 static inline float angle_err_ri(float r1,float i1,float r2,float i2) {
00877 if ((r1==0 && i1==0) || (r2==0 && i2==0)) return 0;
00878
00879
00880 return fast_acos((r1*r2+i1*i2)/(float)(hypot(r1,i1)*hypot(r2,i2)));
00881 }
00882
00889 static inline int goodf(const float *p_f)
00890 {
00891
00892
00893
00894 #if defined(_WIN32)
00895
00896 if ((((int *) p_f)[0] & 0x7f800000) == 0 ||
00897 (((int *) p_f)[0] & 0x7f800000) == 255) {
00898 return 0;
00899 }
00900 #else
00901 using std::fpclassify;
00902 int i = fpclassify(*p_f);
00903 if ( i == FP_NAN || i == FP_INFINITE) return 0;
00904 #endif //_WIN32 || __APPLE__
00905
00906 return 1;
00907 }
00908
00909 static inline int goodf(const double *p_f)
00910 {
00911
00912
00913
00914 #if defined(_WIN32)
00915
00916 if ((((int *) p_f)[0] & 0x7f800000) == 0 ||
00917 (((int *) p_f)[0] & 0x7f800000) == 255) {
00918 return 0;
00919 }
00920 #else
00921 using std::fpclassify;
00922 int i = fpclassify(*p_f);
00923 if ( i == FP_NAN || i == FP_INFINITE) return 0;
00924 #endif //_WIN32 || __APPLE__
00925
00926 return 1;
00927 }
00928
00929 #ifndef _WIN32
00930 static string recv_broadcast(int port);
00931 #endif //_WIN32
00932
00936 static string get_time_label();
00937
00944 static void set_log_level(int argc, char *argv[]);
00945
00956 static inline float eman_copysign(float a, float b)
00957 {
00958 #ifndef WIN32
00959 return copysign(a, b);
00960 #else
00961 int flip = -1;
00962 if ((a <= 0 && b <= 0) || (a > 0 && b > 0)) {
00963 flip = 1;
00964 }
00965 return a * flip;
00966 #endif
00967 }
00968
00981 static inline float eman_erfc(float x)
00982 {
00983 #ifndef WIN32
00984 return (float)erfc(x);
00985 #else
00986 static double a[] = { -1.26551223, 1.00002368,
00987 0.37409196, 0.09678418,
00988 -0.18628806, 0.27886807,
00989 -1.13520398, 1.48851587,
00990 -0.82215223, 0.17087277
00991 };
00992
00993 double result = 1;
00994 double z = fabs(x);
00995 if (z > 0) {
00996 double t = 1 / (1 + 0.5 * z);
00997 double f1 = t * (a[4] + t * (a[5] + t * (a[6] +
00998 t * (a[7] + t * (a[8] + t * a[9])))));
00999 result = t * exp((-z * z) + a[0] + t * (a[1] + t * (a[2] + t * (a[3] + f1))));
01000
01001 if (x < 0) {
01002 result = 2 - result;
01003 }
01004 }
01005 return (float)result;
01006 #endif
01007 }
01008
01020 static void equation_of_plane(const Vec3f& p1, const Vec3f& p2, const Vec3f& p3, float * plane );
01021
01022
01023
01038 static bool point_is_in_convex_polygon_2d(const Vec2f& p1, const Vec2f& p2, const Vec2f& p3, const Vec2f& p4,const Vec2f& actual_point);
01039
01051 static bool point_is_in_triangle_2d(const Vec2f& p1, const Vec2f& p2, const Vec2f& p3, const Vec2f& actual_point);
01052
01059 static void printMatI3D(MIArray3D& mat,
01060 const string str = string(""),
01061 ostream& out = std::cout);
01068 template<class T> static inline T sgn(T& val) {
01069 return (val > 0) ? T(+1) : T(-1);
01070 }
01071
01072
01073
01074
01075
01076
01077
01078
01079
01080
01081
01082
01083
01095 static float* getBaldwinGridWeights( const int& freq_cutoff, const float& P, const float& r, const float& dfreq = 1, const float& alpha=0.5, const float& beta = 0.2);
01096
01103 static bool IsPower2(int x) {
01104 return ( (x>1) && (x & (x-1))==0 );
01105 }
01106
01107
01108 static void apply_precision(float& value, const float& precision) {
01109 float c = ceilf(value);
01110 float f = (float)fast_floor(value);
01111 if (fabs(value - c) < precision) value = c;
01112 else if (fabs(value - f) < precision) value = f;
01113 }
01114 };
01115 }
01116
01117 #endif
01118
01119