Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members

util.h

Go to the documentation of this file.
00001 
00005 /*
00006  * Author: Steven Ludtke, 04/10/2003 (sludtke@bcm.edu)
00007  * Copyright (c) 2000-2006 Baylor College of Medicine
00008  *
00009  * This software is issued under a joint BSD/GNU license. You may use the
00010  * source code in this file under either license. However, note that the
00011  * complete EMAN2 and SPARX software packages have some GPL dependencies,
00012  * so you are responsible for compliance with the licenses of these packages
00013  * if you opt to use BSD licensing. The warranty disclaimer below holds
00014  * in either instance.
00015  *
00016  * This complete copyright notice must be included in any revised version of the
00017  * source code. Additional authorship citations may be added, but existing
00018  * author citations must be preserved.
00019  *
00020  * This program is free software; you can redistribute it and/or modify
00021  * it under the terms of the GNU General Public License as published by
00022  * the Free Software Foundation; either version 2 of the License, or
00023  * (at your option) any later version.
00024  *
00025  * This program is distributed in the hope that it will be useful,
00026  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00027  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
00028  * GNU General Public License for more details.
00029  *
00030  * You should have received a copy of the GNU General Public License
00031  * along with this program; if not, write to the Free Software
00032  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
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 //#include <boost/math/special_functions/fpclassify.hpp>
00055 #include "vec3.h"
00056 
00057 #ifdef WIN32
00058 #include <windows.h>
00059 #define M_PI 3.14159265358979323846f
00060 //#define MAXPATHLEN (MAX_PATH * 4)
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                 //static void file_unlock(FILE * file);
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 //                      printf("%f\t%f\t%f\n",r1*r2+i1*i2,hypot(r1,i1),hypot(r2,i2));
00879 //                      return acos((r1*r2+i1*i2)/(hypot(r1,i1)*hypot(r2,i2)));
00880                         return fast_acos((r1*r2+i1*i2)/(float)(hypot(r1,i1)*hypot(r2,i2)));             // fast_acos also tolerates values very slightly above 1
00881                 }
00882 
00889                 static inline int goodf(const float *p_f)
00890                 {
00891                         //This is the old way to judge a good float, which cause problem on
00892                         //Fedora Core 64 bit system. Now use isinff() and isnanf() on Linux.
00893                         //#if defined(_WIN32) || defined(__APPLE__)
00894                         #if defined(_WIN32)
00895                         // the first is abnormal zero the second is +-inf or NaN
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                         //This is the old way to judge a good float, which cause problem on
00912                         //Fedora Core 64 bit system. Now use isinff() and isnanf() on Linux.
00913                         //#if defined(_WIN32) || defined(__APPLE__)
00914                         #if defined(_WIN32)
00915                         // the first is abnormal zero the second is +-inf or NaN
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 //              /** Get the isosurface value for 3D image.
01073 //               *
01074 //               * @param[in] image 3D image data
01075 //               * @param[in] surface_value threshhold for isosurface valuse
01076 //               * @param[in] smooth boolean to specify whether the smooth value needed
01077 //               *
01078 //               * @return Dict to wrap the points(float array), and triangles(int array)
01079 //               *
01080 //               * @exception ImageDimensionException 3D image only
01081 //               * */
01082 //              static Dict get_isosurface(EMData * image, float surface_value, bool smooth);
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 /* vim: set ts=4 noet nospell: */

Generated on Fri Apr 30 15:38:57 2010 for EMAN2 by  doxygen 1.3.9.1