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
00456 EMData *calc_ccfx( EMData * const with, int y0 = 0, int y1 = -1, bool nosum = false, bool flip = false,bool usez=false);
00457
00458
00467 EMData *make_rotational_footprint(bool unwrap = true);
00468 EMData *make_rotational_footprint_e1(bool unwrap = true);
00469 EMData *make_rotational_footprint_cmc(bool unwrap = true);
00470
00487 EMData *make_footprint(int type=0);
00488
00489
00502 EMData *calc_mutual_correlation(EMData * with, bool tocorner = false, EMData * filter = 0);
00503
00504
00520 EMData *unwrap(int r1 = -1, int r2 = -1, int xs = -1, int dx = 0,
00521 int dy = 0, bool do360 = false, bool weight_radial=true) const;
00522
00523 EMData * unwrap_largerR(int r1,int r2,int xs, float rmax_f);
00524
00525 EMData *oneDfftPolar(int size, float rmax, float MAXR);
00526
00527
00535 void apply_radial_func(float x0, float dx, vector < float >array, bool interp = true);
00536
00537
00549 vector < float >calc_radial_dist(int n, float x0, float dx,bool inten);
00550
00551
00564 vector < float >calc_radial_dist(int n, float x0, float dx, int nwedge, float offset,bool inten);
00565
00566
00570 void cconj();
00571
00572
00580 void add_incoherent(EMData * obj);
00581
00582
00594 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);
00595
00596
00608 vector<float> calc_az_dist(int n, float a0, float da, float rmin,
00609 float rmax);
00610
00611 #if 0
00612 void calc_rcf(EMData * with, vector < float >&sum_array);
00613 #endif
00614
00627 float calc_dist(EMData * second_img, int y_index = 0) const;
00628
00645 EMData *calc_flcf(EMData * with);
00646
00665 EMData *calc_fast_sigma_image( EMData* mask);
00666
00672 EMData *convolute(EMData * with);
00673
00674 #if 0
00675 void create_ctf_map(CtfMapType type, XYData * sf = 0);
00676 #endif
00677
00678
00700 void common_lines(EMData * image1, EMData * image2, int mode = 0,
00701 int steps = 180, bool horizontal = false);
00702
00713 void common_lines_real(EMData * image1, EMData * image2,
00714 int steps = 180, bool horizontal = false);
00715
00729 void cut_slice(const EMData * const map, const Transform& tr, bool interpolate = true);
00730
00745 void uncut_slice(EMData * const map, const Transform& tr) const;
00746
00755 EMData *extract_box(const Transform& cs, const Region& r);
00756
00760 int getResolution() const {
00761 int resolution = 0;
00762 int num = 1;
00763 while(num < get_xsize()) {
00764 resolution++;
00765 num = 1 << resolution;
00766 }
00767
00768 return resolution;
00769 }
00770
00773 void debug_print_parms()
00774 {
00775 std::cout << "Printing EMData params" << std::endl;
00776 for ( Dict::const_iterator it = attr_dict.begin(); it != attr_dict.end(); ++it )
00777 {
00778 std::cout << (it->first) << " " << (it->second).to_str() << std::endl;
00779 }
00780 std::cout << "Done printing EMData params" << std::endl;
00781 }
00782
00788 void set_xyz_origin(float origin_x, float origin_y, float origin_z);
00789
00795 EMData* compute_missingwedge(float wedgeangle, float start = 0.05, float stop = 0.5);
00796
00797 static int totalalloc;
00798 private:
00805 enum EMDataFlags {
00806
00807
00808 EMDATA_BUSY = 1 << 3,
00809 EMDATA_HASCTFF = 1 << 4,
00810 EMDATA_NEEDUPD = 1 << 5,
00811
00812 EMDATA_FLIP = 1 << 7,
00813 EMDATA_PAD = 1 << 8,
00814 EMDATA_FFTODD = 1 << 9,
00815 EMDATA_SHUFFLE = 1 << 10,
00816 EMDATA_FH = 1 << 11,
00817 EMDATA_CPU_NEEDS_UPDATE = 1 << 12,
00818 EMDATA_GPU_NEEDS_UPDATE = 1 << 13,
00819 EMDATA_GPU_RO_NEEDS_UPDATE = 1 << 14
00820 };
00821
00822 void update_stat() const;
00823 void save_byteorder_to_dict(ImageIO * imageio);
00824
00825 private:
00827 mutable Dict attr_dict;
00829 mutable float *rdata;
00831 float *supp;
00832
00835
00836
00838 mutable int flags;
00839
00840 int changecount;
00842 int nx, ny, nz, nxy;
00843 size_t nxyz;
00845 int xoff, yoff, zoff;
00846
00848 Vec3f all_translation;
00849
00850
00851 string path;
00852 int pathnum;
00853
00855 mutable EMData* rot_fp;
00856
00857
00858
00859 class ClipInplaceVariables
00860 {
00861 public:
00862 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) :
00863 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),
00864 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),
00865 prv_z_top(0), prv_z_bottom(0), prv_y_back(0), prv_y_front(0), prv_x_left(0), prv_x_right(0)
00866 {
00867 if ( xtrans > 0 ) x_iter -= xtrans;
00868 if ( x_iter < 0 ) x_iter = 0;
00869 if ( ytrans > 0 ) y_iter -= ytrans;
00870 if ( y_iter < 0 ) y_iter = 0;
00871 if ( ztrans > 0 ) z_iter -= ztrans;
00872 if ( z_iter < 0 ) z_iter = 0;
00873
00874
00875
00876
00877 if ( (new_nz + ztrans) > prv_nz ) new_z_top = new_nz + ztrans - prv_nz;
00878 if ( (new_ny + ytrans) > prv_ny ) new_y_back = new_ny + ytrans - prv_ny;
00879 if ( (new_nx + xtrans) > prv_nx ) new_x_right = new_nx + xtrans - prv_nx;
00880
00881 if ( (new_nz + ztrans) < prv_nz )
00882 {
00883 prv_z_top = prv_nz - new_nz - ztrans;
00884 z_iter -= prv_z_top;
00885 }
00886 if ( (new_ny + ytrans) < prv_ny )
00887 {
00888 prv_y_back = prv_ny - new_ny - ytrans;
00889 y_iter -= prv_y_back;
00890 }
00891 if ( (new_nx + xtrans) < prv_nx )
00892 {
00893 prv_x_right = prv_nx - new_nx - xtrans;
00894 x_iter -= prv_x_right;
00895 }
00896
00897 if ( xtrans > 0 ) prv_x_left = xtrans;
00898 if ( ytrans > 0 ) prv_y_front = ytrans;
00899 if ( ztrans > 0 ) prv_z_bottom = ztrans;
00900
00901 if ( xtrans < 0 ) new_x_left = -xtrans;
00902 if ( ytrans < 0 ) new_y_front = -ytrans;
00903 if ( ztrans < 0 ) new_z_bottom = -ztrans;
00904
00905 }
00906 ~ClipInplaceVariables() {}
00907
00908 int prv_nx, prv_ny, prv_nz, new_nx, new_ny, new_nz;
00909 int xshift, yshift, zshift;
00910 int x_iter, y_iter, z_iter;
00911 int new_z_top, new_z_bottom, new_y_back, new_y_front, new_x_left, new_x_right;
00912 int prv_z_top, prv_z_bottom, prv_y_back, prv_y_front, prv_x_left, prv_x_right;
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
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 EMData * operator/(const EMData & em, float n);
00991
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 EMData * operator/(float n, const EMData & em);
00996
00997 EMData * rsub(const EMData & em, float n);
00998 EMData * rdiv(const EMData & em, float n);
00999
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 EMData * operator/(const EMData & a, const EMData & b);
01004
01005
01006
01007
01008
01009
01010
01011
01012
01013
01014
01015
01016 }
01017
01018
01019 #endif
01020
01021