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__emdata_h__
00037 #define eman__emdata_h__ 1
00038
00039 #ifdef _WIN32
00040 #pragma warning(disable:4819)
00041 #endif //_WIN32
00042
00043 #include <cfloat>
00044 #include <complex>
00045 #include <fstream>
00046
00047 #include "sparx/fundamentals.h"
00048 #include "emutil.h"
00049 #include "util.h"
00050 #include "sparx/emarray.h"
00051 #include "geometry.h"
00052 #include "transform.h"
00053 #ifdef EMAN2_USING_CUDA
00054 #include <cuda_runtime_api.h>
00055 #include "cuda/cuda_util.h"
00056 #endif // EMAN2_USING_CUDA
00057 using std::string;
00058 using std::vector;
00059 using std::map;
00060
00061 using std::ostream;
00062
00063 #include <utility>
00064 using std::pair;
00065
00066 namespace EMAN
00067 {
00068 class ImageIO;
00069 class Ctf;
00070 class XYData;
00071 class Transform;
00072 class GLUtil;
00073
00074 typedef boost::multi_array_ref<float, 2> MArray2D;
00075 typedef boost::multi_array_ref<float, 3> MArray3D;
00076 typedef boost::multi_array_ref<std::complex<float>, 2> MCArray2D;
00077 typedef boost::multi_array_ref<std::complex<float>, 3> MCArray3D;
00078 typedef boost::multi_array<int, 2> MIArray2D;
00079 typedef boost::multi_array<int, 3> MIArray3D;
00080
00087 class EMData
00088 {
00089 friend class GLUtil;
00090
00092 #include "emdata_io.h"
00093
00095 #include "emdata_metadata.h"
00096
00098 #include "emdata_modular.h"
00099
00101 #include "emdata_transform.h"
00102
00104 #include "emdata_core.h"
00105
00107 #include "sparx/emdata_sparx.h"
00108 #ifdef EMAN2_USING_CUDA
00109
00110 #include "emdata_cuda.h"
00111 #endif // EMAN2_USING_CUDA
00112
00113 public:
00114 enum FFTPLACE { FFT_OUT_OF_PLACE, FFT_IN_PLACE };
00115 enum WINDOWPLACE { WINDOW_OUT_OF_PLACE, WINDOW_IN_PLACE };
00116
00118 EMData();
00119 ~ EMData();
00120
00124 explicit EMData(const string& filename, int image_index=0);
00125
00133 EMData(int nx, int ny, int nz=1, bool is_real=true);
00134
00144 EMData(float* data, const int nx, const int ny, const int nz, const Dict& attr_dict = Dict());
00145
00156 EMData(float* data, float* cudadata, const int nx, const int ny, const int nz, const Dict& attr_dict = Dict());
00157
00158
00164 EMData(const EMData& that);
00165
00171 EMData& operator=(const EMData& that);
00172
00173
00181 EMData *get_clip(const Region & area, const float fill = 0) const;
00182
00189 void clip_inplace(const Region & area,const float& fill_value=0);
00190
00195 EMData *get_top_half() const;
00196
00197
00206 EMData *get_rotated_clip(const Transform & xform, const IntSize &size, float scale=1.0);
00207
00218 EMData* window_center(int l);
00219
00220
00235 float *setup4slice(bool redo = true);
00236
00237
00241 void scale(float scale_factor);
00242
00243
00249 void translate(float dx, float dy, float dz);
00250
00251
00255 void translate(const Vec3f &translation);
00256
00257
00264 void translate(int dx, int dy, int dz);
00265
00266
00271 void translate(const Vec3i &translation);
00272
00273
00278 void rotate(const Transform & t);
00279
00280 float max_3D_pixel_error(const Transform &t1, const Transform &t2, float r);
00281
00288 void rotate(float az, float alt, float phi);
00289
00290
00295
00296
00300 inline void transform(const Transform& t) {
00301 ENTERFUNC;
00302 process_inplace("xform",Dict("transform",(Transform*)(&t)));
00303
00304 EXITFUNC;
00305 }
00306
00311 inline void rotate_translate(const Transform & t) {
00312 cout << "Deprecation warning. Please consider using EMData::transform() instead " << endl;
00313 transform(t); }
00314
00324 void rotate_translate(float az, float alt, float phi, float dx, float dy, float dz);
00325
00326
00339 void rotate_translate(float az, float alt, float phi, float dx, float dy,
00340 float dz, float pdx, float pdy, float pdz);
00341
00342
00348 void rotate_x(int dx);
00349
00350
00355 inline void rotate_180() {
00356 ENTERFUNC;
00357 process_inplace("math.rotate.180",Dict());
00358 EXITFUNC;
00359 }
00360
00361
00375 double dot_rotate_translate(EMData * with, float dx, float dy, float da,const bool mirror=false);
00376
00377
00388 EMData *little_big_dot(EMData * little_img, bool do_sigma = false);
00389
00390
00403 EMData *do_radon();
00404
00405
00421 EMData *calc_ccf(EMData * with = 0, fp_flag fpflag = CIRCULANT, bool center=false);
00422
00435 void zero_corner_circulant(const int radius = 0);
00436
00455 EMData *calc_ccfx( EMData * const with, int y0 = 0, int y1 = -1, bool nosum = false, bool flip = false);
00456
00457
00466 EMData *make_rotational_footprint(bool unwrap = true);
00467 EMData *make_rotational_footprint_e1(bool unwrap = true);
00468 EMData *make_rotational_footprint_cmc(bool unwrap = true);
00469
00486 EMData *make_footprint(int type=0);
00487
00488
00501 EMData *calc_mutual_correlation(EMData * with, bool tocorner = false, EMData * filter = 0);
00502
00503
00519 EMData *unwrap(int r1 = -1, int r2 = -1, int xs = -1, int dx = 0,
00520 int dy = 0, bool do360 = false, bool weight_radial=true) const;
00521
00522 EMData * unwrap_largerR(int r1,int r2,int xs, float rmax_f);
00523
00524 EMData *oneDfftPolar(int size, float rmax, float MAXR);
00525
00526
00534 void apply_radial_func(float x0, float dx, vector < float >array, bool interp = true);
00535
00536
00548 vector < float >calc_radial_dist(int n, float x0, float dx,bool inten);
00549
00550
00563 vector < float >calc_radial_dist(int n, float x0, float dx, int nwedge, float offset,bool inten);
00564
00565
00569 void cconj();
00570
00571
00579 void add_incoherent(EMData * obj);
00580
00581
00593 vector <float> calc_hist(int hist_size = 128, float hist_min = 0, float hist_max = 0, const float& brt = 0.0f, const float& cont = 1.0f);
00594
00595
00607 vector<float> calc_az_dist(int n, float a0, float da, float rmin,
00608 float rmax);
00609
00610 #if 0
00611 void calc_rcf(EMData * with, vector < float >&sum_array);
00612 #endif
00613
00626 float calc_dist(EMData * second_img, int y_index = 0) const;
00627
00644 EMData *calc_flcf(EMData * with);
00645
00664 EMData *calc_fast_sigma_image( EMData* mask);
00665
00671 EMData *convolute(EMData * with);
00672
00673 #if 0
00674 void create_ctf_map(CtfMapType type, XYData * sf = 0);
00675 #endif
00676
00677
00699 void common_lines(EMData * image1, EMData * image2, int mode = 0,
00700 int steps = 180, bool horizontal = false);
00701
00712 void common_lines_real(EMData * image1, EMData * image2,
00713 int steps = 180, bool horizontal = false);
00714
00728 void cut_slice(const EMData * const map, const Transform& tr, bool interpolate = true);
00729
00744 void uncut_slice(EMData * const map, const Transform& tr) const;
00745
00754 EMData *extract_box(const Transform& cs, const Region& r);
00755
00759 int getResolution() const {
00760 int resolution = 0;
00761 int num = 1;
00762 while(num < get_xsize()) {
00763 resolution++;
00764 num = 1 << resolution;
00765 }
00766
00767 return resolution;
00768 }
00769
00772 void debug_print_parms()
00773 {
00774 std::cout << "Printing EMData params" << std::endl;
00775 for ( Dict::const_iterator it = attr_dict.begin(); it != attr_dict.end(); ++it )
00776 {
00777 std::cout << (it->first) << " " << (it->second).to_str() << std::endl;
00778 }
00779 std::cout << "Done printing EMData params" << std::endl;
00780 }
00781
00787 void set_xyz_origin(float origin_x, float origin_y, float origin_z);
00788
00794 EMData* compute_missingwedge(float wedgeangle, float start = 0.05, float stop = 0.5);
00795
00796 static int totalalloc;
00797 private:
00804 enum EMDataFlags {
00805
00806
00807 EMDATA_BUSY = 1 << 3,
00808 EMDATA_HASCTFF = 1 << 4,
00809 EMDATA_NEEDUPD = 1 << 5,
00810
00811 EMDATA_FLIP = 1 << 7,
00812 EMDATA_PAD = 1 << 8,
00813 EMDATA_FFTODD = 1 << 9,
00814 EMDATA_SHUFFLE = 1 << 10,
00815 EMDATA_FH = 1 << 11,
00816 EMDATA_CPU_NEEDS_UPDATE = 1 << 12,
00817 EMDATA_GPU_NEEDS_UPDATE = 1 << 13,
00818 EMDATA_GPU_RO_NEEDS_UPDATE = 1 << 14
00819 };
00820
00821 void update_stat() const;
00822 void save_byteorder_to_dict(ImageIO * imageio);
00823
00824 private:
00826 mutable Dict attr_dict;
00828 mutable float *rdata;
00830 float *supp;
00831
00834
00835
00837 mutable int flags;
00838
00839 int changecount;
00841 int nx, ny, nz, nxy;
00842 size_t nxyz;
00844 int xoff, yoff, zoff;
00845
00847 Vec3f all_translation;
00848
00849
00850 string path;
00851 int pathnum;
00852
00854 mutable EMData* rot_fp;
00855
00856
00857
00858 class ClipInplaceVariables
00859 {
00860 public:
00861 ClipInplaceVariables(const int p_nx, const int p_ny, const int p_nz, const int n_nx, const int n_ny, const int n_nz,const int xtrans, const int ytrans, const int ztrans) :
00862 prv_nx(p_nx), prv_ny(p_ny), prv_nz(p_nz), new_nx(n_nx), new_ny(n_ny), new_nz(n_nz), xshift(xtrans), yshift(ytrans), zshift(ztrans),
00863 x_iter(prv_nx), y_iter(prv_ny), z_iter(prv_nz), new_z_top(0), new_z_bottom(0), new_y_back(0), new_y_front(0),new_x_left(0), new_x_right(0),
00864 prv_z_top(0), prv_z_bottom(0), prv_y_back(0), prv_y_front(0), prv_x_left(0), prv_x_right(0)
00865 {
00866 if ( xtrans > 0 ) x_iter -= xtrans;
00867 if ( x_iter < 0 ) x_iter = 0;
00868 if ( ytrans > 0 ) y_iter -= ytrans;
00869 if ( y_iter < 0 ) y_iter = 0;
00870 if ( ztrans > 0 ) z_iter -= ztrans;
00871 if ( z_iter < 0 ) z_iter = 0;
00872
00873
00874
00875
00876 if ( (new_nz + ztrans) > prv_nz ) new_z_top = new_nz + ztrans - prv_nz;
00877 if ( (new_ny + ytrans) > prv_ny ) new_y_back = new_ny + ytrans - prv_ny;
00878 if ( (new_nx + xtrans) > prv_nx ) new_x_right = new_nx + xtrans - prv_nx;
00879
00880 if ( (new_nz + ztrans) < prv_nz )
00881 {
00882 prv_z_top = prv_nz - new_nz - ztrans;
00883 z_iter -= prv_z_top;
00884 }
00885 if ( (new_ny + ytrans) < prv_ny )
00886 {
00887 prv_y_back = prv_ny - new_ny - ytrans;
00888 y_iter -= prv_y_back;
00889 }
00890 if ( (new_nx + xtrans) < prv_nx )
00891 {
00892 prv_x_right = prv_nx - new_nx - xtrans;
00893 x_iter -= prv_x_right;
00894 }
00895
00896 if ( xtrans > 0 ) prv_x_left = xtrans;
00897 if ( ytrans > 0 ) prv_y_front = ytrans;
00898 if ( ztrans > 0 ) prv_z_bottom = ztrans;
00899
00900 if ( xtrans < 0 ) new_x_left = -xtrans;
00901 if ( ytrans < 0 ) new_y_front = -ytrans;
00902 if ( ztrans < 0 ) new_z_bottom = -ztrans;
00903
00904 }
00905 ~ClipInplaceVariables() {}
00906
00907 int prv_nx, prv_ny, prv_nz, new_nx, new_ny, new_nz;
00908 int xshift, yshift, zshift;
00909 int x_iter, y_iter, z_iter;
00910 int new_z_top, new_z_bottom, new_y_back, new_y_front, new_x_left, new_x_right;
00911 int prv_z_top, prv_z_bottom, prv_y_back, prv_y_front, prv_x_left, prv_x_right;
00912 };
00913
00914
00915
00916
00917
00918
00919
00920
00921
00922
00923
00924
00925
00926
00927
00928
00929
00930
00931
00932
00933
00934
00935
00936
00937
00938
00939
00940
00941
00942
00943
00944
00945
00946
00947
00948
00949
00950
00951
00952
00953
00954
00955
00956
00957
00958
00959
00960
00961
00962
00963
00964
00965
00966
00967
00968
00969
00970
00971
00972
00973
00974
00975
00976
00977
00978
00979
00980
00981
00982
00983 };
00984
00985
00986 EMData * operator+(const EMData & em, float n);
00987 EMData * operator-(const EMData & em, float n);
00988 EMData * operator*(const EMData & em, float n);
00989 EMData * operator/(const EMData & em, float n);
00990
00991 EMData * operator+(float n, const EMData & em);
00992 EMData * operator-(float n, const EMData & em);
00993 EMData * operator*(float n, const EMData & em);
00994 EMData * operator/(float n, const EMData & em);
00995
00996 EMData * rsub(const EMData & em, float n);
00997 EMData * rdiv(const EMData & em, float n);
00998
00999 EMData * operator+(const EMData & a, const EMData & b);
01000 EMData * operator-(const EMData & a, const EMData & b);
01001 EMData * operator*(const EMData & a, const EMData & b);
01002 EMData * operator/(const EMData & a, const EMData & b);
01003
01004
01005
01006
01007
01008
01009
01010
01011
01012
01013
01014
01015 }
01016
01017
01018 #endif
01019
01020