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 #include <iomanip>
00037
00038 #include "emdata.h"
00039 #include "all_imageio.h"
00040 #include "ctf.h"
00041
00042 #include <iostream>
00043 using std::cout;
00044 using std::endl;
00045
00046 #include <boost/shared_ptr.hpp>
00047 using boost::shared_ptr;
00048
00049 using namespace EMAN;
00050
00051 void EMData::read_image(const string & filename, int img_index, bool nodata,
00052 const Region * region, bool is_3d)
00053 {
00054 ENTERFUNC;
00055
00056 ImageIO *imageio = EMUtil::get_imageio(filename, ImageIO::READ_ONLY);
00057
00058 if (!imageio) {
00059 throw ImageFormatException("cannot create an image io");
00060 }
00061 else {
00062 int err = imageio->read_header(attr_dict, img_index, region, is_3d);
00063 if (err) {
00064 throw ImageReadException(filename, "imageio read header failed");
00065 }
00066 else {
00067 LstIO * myLstIO = dynamic_cast<LstIO *>(imageio);
00068 if(!myLstIO) attr_dict["source_path"] = filename;
00069 attr_dict["source_n"] = img_index;
00070 if (imageio->is_complex_mode()) {
00071 set_complex(true);
00072 set_fftpad(true);
00073 }
00074 if (attr_dict.has_key("is_fftodd") && (int)attr_dict["is_fftodd"] == 1) {
00075 set_fftodd(true);
00076 }
00077 if ((int) attr_dict["is_complex_ri"] == 1) {
00078 set_ri(true);
00079 }
00080 save_byteorder_to_dict(imageio);
00081
00082 nx = attr_dict["nx"];
00083 ny = attr_dict["ny"];
00084 nz = attr_dict["nz"];
00085 attr_dict.erase("nx");
00086 attr_dict.erase("ny");
00087 attr_dict.erase("nz");
00088
00089 if (!nodata) {
00090
00091 if (region) {
00092 nx = (int)region->get_width();
00093 if (nx <= 0) nx = 1;
00094 ny = (int)region->get_height();
00095 if (ny <= 0) ny = 1;
00096 nz = (int)region->get_depth();
00097 if (nz <= 0) nz = 1;
00098 set_size(nx,ny,nz);
00099 to_zero();
00100 }
00101 else {
00102 set_size(nx, ny, nz);
00103 }
00104
00105
00106
00107
00108 int err = imageio->read_data(get_data(), img_index, region, is_3d);
00109 if (err) {
00110 throw ImageReadException(filename, "imageio read data failed");
00111 }
00112 else {
00113 update();
00114 }
00115 }
00116 else {
00117 if (rdata!=0) EMUtil::em_free(rdata);
00118 rdata=0;
00119 }
00120
00121 }
00122 }
00123
00124 #ifndef IMAGEIO_CACHE
00125 if( imageio )
00126 {
00127 #ifdef HDFIO_CACHE
00128 if(dynamic_cast<HdfIO2*>(imageio)==NULL && dynamic_cast<HdfIO*>(imageio)==NULL) {
00129 #endif //HDFIO_CACHE
00130 delete imageio;
00131 imageio = 0;
00132 #ifdef HDFIO_CACHE
00133 }
00134 #endif //HDFIO_CACHE
00135 }
00136 #endif //IMAGEIO_CACHE
00137
00138 EXITFUNC;
00139 }
00140
00141 void EMData::read_binedimage(const string & filename, int img_index, int binfactor, bool fast, bool is_3d)
00142 {
00143 ENTERFUNC;
00144
00145 ImageIO *imageio = EMUtil::get_imageio(filename, ImageIO::READ_ONLY);
00146
00147 if (!imageio) {
00148 throw ImageFormatException("cannot create an image io");
00149 }
00150 else {
00151 int err = imageio->read_header(attr_dict, img_index, 0, is_3d);
00152 if (err) {
00153 throw ImageReadException(filename, "imageio read header failed");
00154 }
00155 else {
00156 attr_dict["source_path"] = filename;
00157 attr_dict["source_n"] = img_index;
00158 if (imageio->is_complex_mode()) {
00159 set_complex(true);
00160 set_fftpad(true);
00161 }
00162 if (attr_dict.has_key("is_fftodd") && (int)attr_dict["is_fftodd"] == 1) {
00163 set_fftodd(true);
00164 }
00165 if ((int) attr_dict["is_complex_ri"] == 1) {
00166 set_ri(true);
00167 }
00168 save_byteorder_to_dict(imageio);
00169
00170 int ori_nx = nx = attr_dict["nx"];
00171 int ori_ny = ny = attr_dict["ny"];
00172 int ori_nz = nz = attr_dict["nz"];
00173 attr_dict.erase("nx");
00174 attr_dict.erase("ny");
00175 attr_dict.erase("nz");
00176
00177
00178 set_size(nx/binfactor, ny/binfactor, nz/binfactor);
00179
00180
00181 EMData* tempdata = new EMData();
00182 size_t sizeofslice = nx*ny*sizeof(float);
00183
00184
00185 int zbin = binfactor;
00186 if(fast) zbin = 1;
00187
00188 float percent = 0.1f;
00189 for(int k = 0; k < ori_nz; k+=binfactor){
00190 if(k > ori_nz*percent){
00191 cout << float(k)/float(ori_nz) << "% Done!" << endl;
00192 percent+=0.1f;
00193 }
00194
00195 const Region* binregion = new Region(0,0,k,ori_nx,ori_ny,zbin);
00196 tempdata->read_image(filename, 0, false, binregion);
00197
00198 if (binfactor > 1) tempdata->process_inplace("math.meanshrink",Dict("n",binfactor));
00199 size_t offset = nx*ny*k/binfactor;
00200
00201 EMUtil::em_memcpy(get_data()+offset,tempdata->get_data(),sizeofslice);
00202 delete binregion;
00203 }
00204
00205 delete tempdata;
00206 update();
00207 }
00208 }
00209
00210 #ifndef IMAGEIO_CACHE
00211 if( imageio )
00212 {
00213 #ifdef HDFIO_CACHE
00214 if(dynamic_cast<HdfIO2*>(imageio)==NULL && dynamic_cast<HdfIO*>(imageio)==NULL) {
00215 #endif //HDFIO_CACHE
00216 delete imageio;
00217 imageio = 0;
00218 #ifdef HDFIO_CACHE
00219 }
00220 #endif //HDFIO_CACHE
00221 }
00222 #endif //IMAGEIO_CACHE
00223
00224 EXITFUNC;
00225 }
00226
00227 #include <sys/stat.h>
00228
00229 void EMData::write_image(const string & filename, int img_index,
00230 EMUtil::ImageType imgtype,
00231 bool header_only, const Region * region,
00232 EMUtil::EMDataType filestoragetype,
00233 bool use_host_endian)
00234 {
00235 ENTERFUNC;
00236
00237 struct stat fileinfo;
00238 if ( region && stat(filename.c_str(),&fileinfo) != 0 ) throw UnexpectedBehaviorException("To write an image using a region the file must already exist and be the correct dimensions");
00239
00240 if (is_complex() && is_shuffled())
00241 fft_shuffle();
00242
00243 if (imgtype == EMUtil::IMAGE_UNKNOWN) {
00244 const char *ext = strrchr(filename.c_str(), '.');
00245 if (ext) {
00246 ext++;
00247 imgtype = EMUtil::get_image_ext_type(ext);
00248 }
00249 }
00250 ImageIO::IOMode rwmode = ImageIO::READ_WRITE;
00251
00252
00253 attr_dict["nx"] = nx;
00254 attr_dict["ny"] = ny;
00255 attr_dict["nz"] = nz;
00256 attr_dict["changecount"] = changecount;
00257
00258 #ifndef HDFIO_CACHE
00259 if (Util::is_file_exist(filename)) {
00260 LOGVAR("file exists");
00261 if (!header_only && region == 0) {
00262 ImageIO * tmp_imageio = EMUtil::get_imageio(filename, ImageIO::READ_ONLY,
00263 imgtype);
00264 if (tmp_imageio->is_single_image_format()) {
00265 rwmode = ImageIO::WRITE_ONLY;
00266 }
00267
00268 #ifndef IMAGEIO_CACHE
00269 if( tmp_imageio )
00270 {
00271 delete tmp_imageio;
00272 tmp_imageio = 0;
00273 }
00274 #endif //IMAGEIO_CACHE
00275 }
00276 }
00277 #endif //HDFIO_CACHE
00278
00279 LOGVAR("getimageio %d",rwmode);
00280 ImageIO *imageio = EMUtil::get_imageio(filename, rwmode, imgtype);
00281 if (!imageio) {
00282 throw ImageFormatException("cannot create an image io");
00283 }
00284 else {
00285 update_stat();
00286
00287
00288
00289
00290 LOGVAR("header write %d",img_index);
00291
00292 switch(filestoragetype) {
00293 case EMUtil::EM_UINT:
00294 attr_dict["datatype"] = (int)EMUtil::EM_UINT;
00295 break;
00296 case EMUtil::EM_USHORT:
00297 attr_dict["datatype"] = (int)EMUtil::EM_USHORT;
00298 break;
00299 case EMUtil::EM_SHORT:
00300 attr_dict["datatype"] = (int)EMUtil::EM_SHORT;
00301 break;
00302 case EMUtil::EM_CHAR:
00303 attr_dict["datatype"] = (int)EMUtil::EM_CHAR;
00304 break;
00305 case EMUtil::EM_UCHAR:
00306 attr_dict["datatype"] = (int)EMUtil::EM_UCHAR;
00307 break;
00308 default:
00309 attr_dict["datatype"] = (int)EMUtil::EM_FLOAT;;
00310 }
00311
00312 int err = imageio->write_header(attr_dict, img_index, region, filestoragetype,
00313 use_host_endian);
00314 if (err) {
00315 throw ImageWriteException(filename, "imageio write header failed");
00316 }
00317 else {
00318 if (!header_only) {
00319 if (imgtype == EMUtil::IMAGE_LST) {
00320 const char *reffile = attr_dict["LST.reffile"];
00321 if (strcmp(reffile, "") == 0) {
00322 reffile = path.c_str();
00323 }
00324 int refn = attr_dict["LST.refn"];
00325 if (refn < 0) {
00326 refn = pathnum;
00327 }
00328
00329 const char *comment = attr_dict["LST.comment"];
00330 char *lstdata = new char[1024];
00331 sprintf(lstdata, "%d\t%s", refn, reffile);
00332 if(strcmp(comment, "") != 0) {
00333 sprintf(lstdata+strlen(lstdata), "\t%s\n", comment);
00334 }
00335 else {
00336 strcat(lstdata, "\n");
00337 }
00338 err = imageio->write_data((float*)lstdata, img_index,
00339 region, filestoragetype, use_host_endian);
00340 if( lstdata )
00341 {
00342 delete [] lstdata;
00343 lstdata = 0;
00344 }
00345 }
00346 if (imgtype == EMUtil::IMAGE_LSTFAST) {
00347 const char *reffile = attr_dict["LST.reffile"];
00348 if (strcmp(reffile, "") == 0) {
00349 reffile = path.c_str();
00350 }
00351 int refn = attr_dict["LST.refn"];
00352 if (refn < 0) {
00353 refn = pathnum;
00354 }
00355
00356 const char *comment = attr_dict["LST.comment"];
00357 char *lstdata = new char[1024];
00358 sprintf(lstdata, "%d\t%s", refn, reffile);
00359 if(strcmp(comment, "") != 0) {
00360 sprintf(lstdata+strlen(lstdata), "\t%s\n", comment);
00361 }
00362 else {
00363 strcat(lstdata, "\n");
00364 }
00365 err = imageio->write_data((float*)lstdata, img_index,
00366 region, filestoragetype, use_host_endian);
00367 if( lstdata )
00368 {
00369 delete [] lstdata;
00370 lstdata = 0;
00371 }
00372 }
00373 else {
00374 err = imageio->write_data(get_data(), img_index, region, filestoragetype,
00375 use_host_endian);
00376 }
00377 if (err) {
00378 imageio->flush();
00379 throw ImageWriteException(filename, "imageio write data failed");
00380 }
00381 }
00382 }
00383 }
00384
00385 if (!(imgtype == EMUtil::IMAGE_PNG)) {
00386 imageio->flush();
00387 }
00388
00389 #ifndef IMAGEIO_CACHE
00390 if( imageio )
00391 {
00392 #ifdef HDFIO_CACHE
00393 if(dynamic_cast<HdfIO2*>(imageio)==NULL && dynamic_cast<HdfIO*>(imageio)==NULL) {
00394 #endif //HDFIO_CACHE
00395 delete imageio;
00396 imageio = 0;
00397 #ifdef HDFIO_CACHE
00398 }
00399 #endif //HDFIO_CACHE
00400 }
00401 #endif //IMAGEIO_CACHE
00402
00403 EXITFUNC;
00404 }
00405
00406
00407 void EMData::append_image(const string & filename,
00408 EMUtil::ImageType imgtype, bool header_only)
00409 {
00410 ENTERFUNC;
00411 write_image(filename, -1, imgtype, header_only, 0);
00412 EXITFUNC;
00413 }
00414
00415
00416 void EMData::write_lst(const string & filename, const string & reffile,
00417 int refn, const string & comment)
00418 {
00419 ENTERFUNC;
00420 attr_dict["LST.reffile"] = reffile;
00421 attr_dict["LST.refn"] = refn;
00422 attr_dict["LST.comment"] = comment;
00423 write_image(filename, -1, EMUtil::IMAGE_LST, false);
00424 EXITFUNC;
00425 }
00426
00427
00428 void EMData::print_image(const string str, ostream& out) {
00429 out << "Printing EMData object: " << str << std::endl;
00430 int nx = get_xsize();
00431 int ny = get_ysize();
00432 int nz = get_zsize();
00433 for (int iz = 0; iz < nz; iz++) {
00434 out << "(z = " << iz << " slice)" << std::endl;
00435 for (int ix = 0; ix < nx; ix++) {
00436 for (int iy = 0; iy < ny; iy++) {
00437 out << setiosflags(std::ios::fixed)
00438 << setiosflags(std::ios_base::scientific)
00439 << std::setw(12)
00440 << std::setprecision(5) << (*this)(ix,iy,iz) << " ";
00441 if (((iy+1) % 6) == 0) {
00442 out << std::endl << " ";
00443 }
00444 }
00445 out << std::endl;
00446 }
00447 }
00448 }
00449
00450 vector < shared_ptr<EMData> > EMData::read_images(const string & filename, vector < int >img_indices,
00451 bool header_only)
00452 {
00453 ENTERFUNC;
00454
00455 int total_img = EMUtil::get_image_count(filename);
00456 size_t num_img = img_indices.size();
00457
00458 for (size_t i = 0; i < num_img; i++) {
00459 if (img_indices[i] < 0 && img_indices[i] >= total_img) {
00460 throw OutofRangeException(0, total_img, img_indices[i], "image index");
00461 }
00462 }
00463
00464 size_t n = (num_img == 0 ? total_img : num_img);
00465
00466 vector< shared_ptr<EMData> > v;
00467 for (size_t j = 0; j < n; j++) {
00468 shared_ptr<EMData> d(new EMData());
00469 size_t k = (num_img == 0 ? j : img_indices[j]);
00470 try {
00471 d->read_image(filename, (int)k, header_only);
00472 }
00473 catch(E2Exception &e) {
00474 throw(e);
00475 }
00476 if ( d != 0 )
00477 {
00478 v.push_back(d);
00479 }
00480 else
00481 throw ImageReadException(filename, "imageio read data failed");
00482 }
00483
00484 EXITFUNC;
00485 return v;
00486 }
00487
00488
00489 vector < shared_ptr<EMData> >EMData::read_images_ext(const string & filename, int img_index_start,
00490 int img_index_end, bool header_only,
00491 const string & ext)
00492 {
00493 ENTERFUNC;
00494
00495 if (img_index_end < img_index_start) {
00496 throw InvalidValueException(img_index_end, "image index end < image index start");
00497 }
00498 string new_filename = filename;
00499 new_filename = new_filename.insert(new_filename.rfind("."), ext);
00500 int num_img = EMUtil::get_image_count(new_filename);
00501
00502 if (img_index_start < 0 || img_index_start >= num_img) {
00503 throw OutofRangeException(0, num_img-1, img_index_start, "image index start");
00504 }
00505
00506 if (img_index_end >= num_img) {
00507 img_index_end = num_img - 1;
00508 }
00509
00510 vector < shared_ptr<EMData> >v;
00511
00512 for (int i = img_index_start; i < img_index_end; i++) {
00513 shared_ptr<EMData> d(new EMData());
00514 try {
00515 d->read_image(new_filename, i, header_only);
00516 }
00517 catch(E2Exception &e) {
00518 throw(e);
00519 }
00520 v.push_back(d);
00521 }
00522 EXITFUNC;
00523 return v;
00524 }
00525