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

hdfio2.cpp

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 #ifdef EM_HDF5
00037 
00038 //#define DEBUGHDF      1
00039 
00040 #include "hdfio2.h"
00041 #include "geometry.h"
00042 #include "ctf.h"
00043 #include "emassert.h"
00044 #include "transform.h"
00045 #include "ctf.h"
00046 #include <iostream>
00047 #include <cstring>
00048 
00049 #ifndef WIN32
00050         #include <sys/param.h>
00051 #else
00052         #define  MAXPATHLEN (MAX_PATH * 4)
00053 #endif  //WIN32
00054 
00055 using namespace EMAN;
00056 
00057 static const int ATTR_NAME_LEN = 128;
00058 
00059 HdfIO2::HdfIO2(const string & hdf_filename, IOMode rw)
00060 :       nx(1), ny(1), nz(1), is_exist(false),
00061         file(-1), group(-1), filename(hdf_filename),
00062         rw_mode(rw), initialized(false), rendermin(0), rendermax(0)
00063 {
00064         accprop=H5Pcreate(H5P_FILE_ACCESS);
00065 
00066         //STDIO file driver has 2G size limit on 32 bit Linux system
00067         H5Pset_fapl_sec2( accprop );
00068         //H5Pset_fapl_stdio( accprop );
00069 
00070 //      H5Pset_fapl_core( accprop, 1048576, 0  );
00071 //      H5Pset_cache(accprop)
00072         hsize_t dims=1;
00073         simple_space=H5Screate_simple(1,&dims,NULL);
00074 
00075         meta_attr_dict = Dict();
00076 }
00077 
00078 HdfIO2::~HdfIO2()
00079 {
00080         H5Sclose(simple_space);
00081         H5Pclose(accprop);
00082     if (group >= 0) {
00083         H5Gclose(group);
00084     }
00085     if (file >= 0) {
00086                 H5Fflush(file,H5F_SCOPE_GLOBAL);        // If there were no resource leaks, this wouldn't be necessary...
00087                 H5Fclose(file);
00088     }
00089 #ifdef DEBUGHDF
00090         printf("HDf close\n");
00091 #endif
00092 }
00093 
00094 // This reads an already opened attribute and returns the results as an EMObject
00095 // The attribute is not closed
00096 EMObject HdfIO2::read_attr(hid_t attr) {
00097         hid_t type = H5Aget_type(attr);
00098         hid_t spc = H5Aget_space(attr);
00099         H5T_class_t cls = H5Tget_class(type);
00100         size_t sz = H5Tget_size(type);                                          // storage size, arrays handled in the 'space'
00101         hssize_t pts = H5Sget_simple_extent_npoints(spc);       // number of points > 1 if an array of floats or integers
00102 
00103         EMObject ret(0);
00104         char c;
00105         int i;
00106 //      unsigned int ui;
00107         float f,*fa;
00108         int * ia;
00109 //      unsigned int * uia;
00110         double d;
00111         char *s;
00112         vector <float> fv((size_t)pts);
00113         vector <int> iv((size_t)pts);
00114 //      vector <unsigned int> uiv(pts);
00115 
00116         float *matrix;
00117         Transform* t;
00118         Ctf* ctf;
00119 //      int r, c, k=0;
00120 
00121         switch (cls) {
00122         case H5T_INTEGER:
00123                 if(sz==1) {
00124                         H5Aread(attr,H5T_NATIVE_CHAR,&c);
00125                         bool b = false;
00126                         if(c=='T') {
00127                                 b = true;
00128                         }
00129                         else if(c=='F') {
00130                                 b = false;
00131                         }
00132                         ret = EMObject(b);
00133                 }
00134                 else if(sz==4) {
00135                         if(pts==1) {
00136                                 H5Aread(attr,H5T_NATIVE_INT,&i);
00137                                 ret=EMObject(i);
00138                         }
00139                         else {
00140                                 ia=(int *)malloc((size_t)pts*sizeof(int));
00141                                 H5Aread(attr,H5T_NATIVE_INT,ia);
00142                                 for (i=0; i<pts; i++) iv[i]=ia[i];
00143                                 free(ia);
00144                                 ret=EMObject(iv);
00145                         }
00146                 }
00147                 break;
00148 //      case H5T_UNSIGNED_INTEGER:
00149 //              if(pts==1) {
00150 //                      H5Aread(attr,H5T_NATIVE_UINT,&ui);
00151 //                      ret=EMObject(ui);
00152 //              }
00153 //              else {
00154 //                      uia=(unsigned int *)malloc(pts*sizeof(unsigned int));
00155 //                      H5Aread(attr,H5T_NATIVE_UINT,uia);
00156 //                      for (i=0; i<pts; i++) uiv[i]=uia[i];
00157 //                      free(uia);
00158 //                      ret=EMObject(uiv);
00159 //              }
00160 //              break;
00161         case H5T_FLOAT:
00162                 if (sz==4) {
00163                         if (pts==1) {
00164                                 H5Aread(attr,H5T_NATIVE_FLOAT,&f);
00165                                 ret=EMObject(f);
00166                         }
00167                         else {
00168                                 fa=(float *)malloc((size_t)pts*sizeof(float));
00169                                 H5Aread(attr,H5T_NATIVE_FLOAT,fa);
00170                                 for (i=0; i<pts; i++) fv[i]=fa[i];
00171                                 free(fa);
00172                                 ret=EMObject(fv);
00173                         }
00174                 }
00175                 else if (sz==8) {
00176                         H5Aread(attr,H5T_NATIVE_DOUBLE,&d);
00177                         ret=EMObject(d);
00178                 }
00179                 break;
00180         case H5T_STRING:
00181                 s=(char *)malloc(sz+1);
00182                 H5Aread(attr,type,s);
00183 //              H5Aread(attr,H5T_NATIVE_CHAR,s);
00184                 if(s[0] == 'O' && isdigit(s[1])) {
00185                         ctf = new EMAN1Ctf();
00186                         try {
00187                                 ctf->from_string(string(s));
00188                                 ret = EMObject(ctf);
00189                         }
00190                         catch(...) {
00191                                 ret=EMObject(s);
00192                         }
00193                         delete ctf;
00194                 }
00195                 else if(s[0] == 'E' && isdigit(s[1])) {
00196                         ctf = new EMAN2Ctf();
00197                         try {
00198                                 ctf->from_string(string(s));
00199                                 ret = EMObject(ctf);
00200                         }
00201                         catch(...) {
00202                                 ret=EMObject(s);
00203                         }
00204 
00205                         delete ctf;
00206                 }
00207                 else {
00208                         ret=EMObject(s);
00209                 }
00210                 free(s);
00211                 break;
00212         case H5T_COMPOUND:
00213                 matrix = (float*)malloc(12*sizeof(float));
00214                 H5Aread(attr, type, matrix);
00215 //              ret.create_transform3d_by_array(trans3d);
00216                 t = new Transform(matrix);
00217                 ret = EMObject(t);
00218                 free(matrix);
00219                 delete t; t=0;
00220 
00221 //              trans3d = (float*)malloc(16*sizeof(float));     //16 float for a Transform3D object
00222 //              H5Aread(attr, type, trans3d);
00223 // //           ret.create_transform3d_by_array(trans3d);
00224 //              trans = new Transform3D();
00225 //              for(r=0; r<4; ++r) {
00226 //                      for(c=0; c<4; ++c) {
00227 //                              trans->set(r, c, trans3d[k]);
00228 //                              ++k;
00229 //                      }
00230 //              }
00231 //              ret = EMObject(trans);
00232 //              free(trans3d);
00233                 break;
00234         default:
00235                 LOGERR("Unhandled HDF5 metadata %d", cls);
00236         }
00237 
00238         H5Sclose(spc);
00239         H5Tclose(type);
00240 
00241         return ret;
00242 }
00243 
00244 // This writes an attribute with specified name to a given open object
00245 // The attribute is opened and closed. returns 0 on success
00246 int HdfIO2::write_attr(hid_t loc,const char *name,EMObject obj) {
00247         hid_t type=0;
00248         hid_t spc=0;
00249         hsize_t dims=1;
00250         vector <float> fv;
00251         vector <int> iv;
00252         switch(obj.get_type())
00253         {
00254         case EMObject::BOOL:
00255                 type=H5Tcopy(H5T_NATIVE_CHAR);
00256                 spc=H5Scopy(simple_space);
00257                 break;
00258         case EMObject::SHORT:
00259         case EMObject::INT:
00260                 type=H5Tcopy(H5T_NATIVE_INT);
00261                 spc=H5Scopy(simple_space);
00262                 break;
00263         case EMObject::UNSIGNEDINT:
00264                 type=H5Tcopy(H5T_NATIVE_UINT);
00265                 spc=H5Scopy(simple_space);
00266                 break;
00267         case EMObject::FLOAT:
00268                 type=H5Tcopy(H5T_NATIVE_FLOAT);
00269                 spc=H5Scopy(simple_space);
00270                 break;
00271         case EMObject::DOUBLE:
00272                 type=H5Tcopy(H5T_NATIVE_DOUBLE);
00273                 spc=H5Scopy(simple_space);
00274                 break;
00275         case EMObject::STRING:
00276         case EMObject::CTF:
00277                 type=H5Tcopy(H5T_C_S1);
00278                 H5Tset_size(type,strlen((const char *)obj)+1);
00279                 spc=H5Screate(H5S_SCALAR);
00280                 break;
00281         case EMObject::FLOATARRAY:
00282                 type=H5Tcopy(H5T_NATIVE_FLOAT);
00283                 fv=obj;
00284                 dims=fv.size();
00285                 spc=H5Screate_simple(1,&dims,NULL);
00286                 break;
00287         case EMObject::INTARRAY:
00288                 type=H5Tcopy(H5T_NATIVE_INT);
00289                 iv=obj;
00290                 dims=iv.size();
00291                 spc=H5Screate_simple(1,&dims,NULL);
00292                 break;
00293         case EMObject::TRANSFORM:
00294                 type = H5Tcreate(H5T_COMPOUND, 12 * sizeof(float)); //Transform is a 3x4 matrix
00295                 H5Tinsert(type, "00", 0, H5T_NATIVE_FLOAT);
00296                 H5Tinsert(type, "01", 1*sizeof(float), H5T_NATIVE_FLOAT);
00297                 H5Tinsert(type, "02", 2*sizeof(float), H5T_NATIVE_FLOAT);
00298                 H5Tinsert(type, "03", 3*sizeof(float), H5T_NATIVE_FLOAT);
00299                 H5Tinsert(type, "10", 4*sizeof(float), H5T_NATIVE_FLOAT);
00300                 H5Tinsert(type, "11", 5*sizeof(float), H5T_NATIVE_FLOAT);
00301                 H5Tinsert(type, "12", 6*sizeof(float), H5T_NATIVE_FLOAT);
00302                 H5Tinsert(type, "13", 7*sizeof(float), H5T_NATIVE_FLOAT);
00303                 H5Tinsert(type, "20", 8*sizeof(float), H5T_NATIVE_FLOAT);
00304                 H5Tinsert(type, "21", 9*sizeof(float), H5T_NATIVE_FLOAT);
00305                 H5Tinsert(type, "22", 10*sizeof(float), H5T_NATIVE_FLOAT);
00306                 H5Tinsert(type, "23", 11*sizeof(float), H5T_NATIVE_FLOAT);
00307                 H5Tpack(type);
00308 
00309                 dims = 1;       //one compound type
00310                 spc = H5Screate_simple(1, &dims, NULL);
00311                 break;
00312         case EMObject::TRANSFORMARRAY:
00313         case EMObject::STRINGARRAY:
00314         case EMObject::EMDATA:
00315         case EMObject::XYDATA:
00316         case EMObject::FLOAT_POINTER:
00317         case EMObject::INT_POINTER:
00318         case EMObject::VOID_POINTER:
00319                 return -1;
00320                 break;
00321         case EMObject::UNKNOWN:
00322                 break;
00323         }
00324 
00325     //we need this delete attribute call here, even we called erase_header()
00326     //at the beginning of write_header(), since the  "imageid_max" need be updated correctly.
00327         if( H5Adelete(loc,name) < 0 ) {
00328 #ifdef DEBUGHDF
00329                 LOGERR("Attribute %s deletion error in write_attr().\n", name);
00330 #endif
00331         }
00332         else {
00333 #ifdef DEBUGHDF
00334                 printf("delete attribute %s successfully in write_attr().\n", name);
00335 #endif
00336         }
00337         hid_t attr = H5Acreate(loc,name,type,spc,H5P_DEFAULT);
00338 
00339         bool b;
00340         char c;
00341         int i;
00342         short si;
00343         float f,*fa;
00344         int * ia;
00345         unsigned int ui;
00346         double d;
00347         const char *s;
00348         Transform * tp;
00349         switch(obj.get_type()) {
00350         case EMObject::BOOL:
00351                 b = (bool)obj;
00352                 if(b) {
00353                         c = 'T';
00354                 } else {
00355                         c = 'F';
00356                 }
00357                 H5Awrite(attr,type,&c);
00358                 break;
00359         case EMObject::SHORT:
00360                 si = (short)obj;
00361                 i = (int)si;
00362                 H5Awrite(attr,type,&i);
00363                 break;
00364         case EMObject::INT:
00365                 i=(int)obj;
00366                 H5Awrite(attr,type,&i);
00367                 break;
00368         case EMObject::UNSIGNEDINT:
00369                 ui=(unsigned int)obj;
00370                 H5Awrite(attr,type,&ui);
00371                 break;
00372         case EMObject::FLOAT:
00373                 f=(float)obj;
00374                 H5Awrite(attr,type,&f);
00375                 break;
00376         case EMObject::DOUBLE:
00377                 d=(double)obj;
00378                 H5Awrite(attr,type,&d);
00379                 break;
00380         case EMObject::STRING:
00381         case EMObject::CTF:
00382                 s=(const char *)obj;
00383                 H5Awrite(attr,type,s);
00384                 break;
00385         case EMObject::FLOATARRAY:
00386                 fa=(float *)malloc(fv.size()*sizeof(float));
00387                 for (ui=0; ui<fv.size(); ui++) fa[ui]=fv[ui];
00388                 H5Awrite(attr,type,fa);
00389                 free(fa);
00390                 break;
00391         case EMObject::INTARRAY:
00392                 ia=(int *)malloc(iv.size()*sizeof(int));
00393                 for (ui=0; ui<iv.size(); ui++) ia[ui]=iv[ui];
00394                 H5Awrite(attr,type,ia);
00395                 free(ia);
00396                 break;
00397         case EMObject::TRANSFORM:
00398         {
00399                 tp = (Transform *)obj;
00400                 fa = (float *)malloc(12*sizeof(float));
00401                 int r, c, k=0;
00402                 for(r=0; r<3; ++r) {
00403                         for(c=0; c<4; ++c) {
00404                                 fa[k] = tp->at(r,c);
00405                                 ++k;
00406                         }
00407                 }
00408                 H5Awrite(attr,type,fa);
00409                 free(fa);
00410         }
00411                 break;
00412 //      case EMObject::STRINGARRAY:
00413 //      case EMObject::EMDATA:
00414 //      case EMObject::XYDATA:
00415 //              return -1;
00416 //              break;
00417         default:
00418                 LOGERR("Unhandled HDF5 metadata '%s'", name);
00419 
00420         }
00421 
00422         herr_t ret1 = H5Tclose(type);
00423         herr_t ret2 = H5Sclose(spc);
00424         herr_t ret3 = H5Aclose(attr);
00425         if(ret1>=0 && ret2>=0 && ret3>=0) {
00426                 return 0;
00427         }
00428         else {
00429                 LOGERR("close error in write_attr()\n");
00430                 return -1;
00431         }
00432 }
00433 
00434 // Initializes the file for read-only or read-write access
00435 // Data is stored under /MDF/images
00436 // An attribute named imageid_max stores the number of the highest
00437 // numbered image in the file.
00438 // A group is then made for each individual image, all metadata for the
00439 // individual images is currently associated with the GROUP, not the dataset
00440 // dataset-specific data could also be associated with the dataset in
00441 // future. At the moment, there is only a single dataset in each group.
00442 void HdfIO2::init()
00443 {
00444         ENTERFUNC;
00445         if (initialized) {
00446                 return;
00447         }
00448 #ifdef DEBUGHDF
00449         printf("init\n");
00450 #endif
00451 
00452         H5Eset_auto(0, 0);      // Turn off console error logging.
00453 
00454         if (rw_mode == READ_ONLY) {
00455                 file = H5Fopen(filename.c_str(), H5F_ACC_RDONLY, accprop);
00456                 if (file<0) throw FileAccessException(filename);
00457         }
00458         else {
00459                 file = H5Fopen(filename.c_str(), H5F_ACC_RDWR, accprop);
00460                 if (file < 0) {
00461                         file = H5Fcreate(filename.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, accprop);
00462                         if (file < 0) {
00463                                 throw FileAccessException(filename);
00464                         }
00465                         else {
00466 #ifdef DEBUGHDF
00467                                 printf("File truncated or new file created\n");
00468 #endif
00469                         }
00470                 }
00471         }
00472 
00473         group=H5Gopen(file,"/MDF/images");
00474         if (group<0) {
00475                 if (rw_mode == READ_ONLY) throw ImageReadException(filename,"HDF5 file has no image data (no /MDF group)");
00476                 group=H5Gcreate(file,"/MDF",64);                // create the group for Macromolecular data
00477                 if (group<0) throw ImageWriteException(filename,"Unable to add image group (/MDF) to HDF5 file");
00478                 H5Gclose(group);
00479                 group=H5Gcreate(file,"/MDF/images",4096);               // create the group for images/volumes
00480                 if (group<0) throw ImageWriteException(filename,"Unable to add image group (/MDF/images) to HDF5 file");
00481                 write_attr(group,"imageid_max",EMObject(-1));
00482         }
00483         else {  //read the meta attributes for all images
00484                 int nattr=H5Aget_num_attrs(group);
00485 
00486                 char name[ATTR_NAME_LEN];
00487                 for (int i=0; i<nattr; i++) {
00488                         hid_t attr=H5Aopen_idx(group, i);
00489                         H5Aget_name(attr,127,name);
00490 
00491                         EMObject val=read_attr(attr);
00492                         meta_attr_dict["DDD."+string(name)]=val;
00493 
00494                         H5Aclose(attr);
00495                 }
00496 
00497         }
00498         initialized = true;
00499         EXITFUNC;
00500 }
00501 
00502 
00503 // If this version of init() returns -1, then we have an old-style HDF5 file
00504 int HdfIO2::init_test()
00505 {
00506         ENTERFUNC;
00507         if (initialized) {
00508                 return 1;
00509         }
00510 #ifdef DEBUGHDF
00511         printf("init_test\n");
00512 #endif
00513 
00514         H5Eset_auto(0, 0);      // Turn off console error logging.
00515 
00516         hid_t fileid = H5Fopen(filename.c_str(), H5F_ACC_RDONLY, H5Pcreate(H5P_FILE_ACCESS));
00517         hid_t groupid = H5Gopen(fileid, "/");
00518         hid_t attid = H5Aopen_name(groupid, "num_dataset");
00519 
00520         if (attid < 0) {
00521                 H5Gclose(groupid);
00522                 H5Fclose(fileid);
00523                 init();
00524                 EXITFUNC;
00525                 return 0;
00526         }
00527         else {
00528                 H5Aclose(attid);
00529                 H5Gclose(groupid);
00530                 H5Fclose(fileid);
00531                 EXITFUNC;
00532                 return -1;
00533         }
00534 }
00535 
00536 bool HdfIO2::is_valid(const void *first_block)
00537 {
00538         ENTERFUNC;
00539 
00540         if (first_block) {
00541                 char signature[8] = { 137,72,68,70,13,10,26,10 };
00542                 if (strncmp((const char *)first_block,signature,8)==0) return true;
00543                 // const char* f=(const char *)first_block;
00544                 // printf("bad hdf signature %d %d %d %d %d %d %d %d",f[0],f[1],f[2],f[3],f[4],f[5],f[6],f[7]);
00545                 return false;
00546         }
00547         EXITFUNC;
00548         return false;
00549 }
00550 
00551 // Reads all of the attributes from the /MDF/images/<imgno> group
00552 int HdfIO2::read_header(Dict & dict, int image_index, const Region * area, bool)
00553 {
00554         ENTERFUNC;
00555         init();
00556 
00558         size_t meta_attr_size = meta_attr_dict.size();
00559         if(meta_attr_size!=0) {
00560                 for (size_t i=0; i<meta_attr_size; ++i) {
00561                         dict[meta_attr_dict.keys()[i]] = meta_attr_dict.values()[i];
00562                 }
00563         }
00564 
00565 #ifdef DEBUGHDF
00566         printf("read_head %d\n", image_index);
00567 #endif
00568         int i;
00569         // Each image is in a group for later expansion. Open the group
00570         char ipath[50];
00571         sprintf(ipath,"/MDF/images/%d", image_index);
00572         hid_t igrp=H5Gopen(file, ipath);
00573 
00574         int nattr=H5Aget_num_attrs(igrp);
00575 
00576         char name[ATTR_NAME_LEN];
00577         for (i=0; i<nattr; i++) {
00578                 hid_t attr=H5Aopen_idx(igrp, i);
00579                 H5Aget_name(attr,127,name);
00580                 if (strncmp(name,"EMAN.", 5)!=0) {
00581                         H5Aclose(attr);
00582                         continue;
00583                 }
00584                 try {
00585                         EMObject val=read_attr(attr);
00586                         dict[name+5]=val;
00587                 }
00588                 catch(...) {
00589                         printf("Error reading HDF attribute %s\n",name+5);
00590                 }
00591                 H5Aclose(attr);
00592         }
00593 
00594         if(dict.has_key("ctf")) {
00595                 string ctfString = (string)dict["ctf"];
00596                 if(ctfString.substr(0, 1) == "O") {
00597                         Ctf * ctf_ = new EMAN1Ctf();
00598                         ctf_->from_string(ctfString);
00599                         dict.erase("ctf");
00600                         dict["ctf"] = ctf_;
00601                         delete ctf_;
00602                 }
00603                 else if(ctfString.substr(0, 1) == "E") {
00604                         Ctf * ctf_ = new EMAN2Ctf();
00605                         ctf_->from_string(ctfString);
00606                         dict.erase("ctf");
00607                         dict["ctf"] = ctf_;
00608                         delete ctf_;
00609                 }
00610         }
00611 
00612         if(area) {
00613                 check_region(area, IntSize(dict["nx"], dict["ny"], dict["nz"]), false, false);
00614 
00615                 dict["nx"] = area->get_width();
00616                 dict["ny"] = area->get_height();
00617                 dict["nz"] = area->get_depth();
00618 
00619                 if( dict.has_key("apix_x") && dict.has_key("apix_y") && dict.has_key("apix_z") )
00620                 {
00621                         if( dict.has_key("origin_x") && dict.has_key("origin_y") && dict.has_key("origin_z") )
00622                         {
00623                                 float xorigin = dict["origin_x"];
00624                                 float yorigin = dict["origin_y"];
00625                                 float zorigin = dict["origin_z"];
00626 
00627                                 float apix_x = dict["apix_x"];
00628                                 float apix_y = dict["apix_y"];
00629                                 float apix_z = dict["apix_z"];
00630 
00631                                 dict["origin_x"] = xorigin + apix_x * area->origin[0];
00632                                 dict["origin_y"] = yorigin + apix_y * area->origin[1];
00633                                 dict["origin_z"] = zorigin + apix_z * area->origin[2];
00634                         }
00635                 }
00636         }
00637 
00638         H5Gclose(igrp);
00639 
00640         //Get the data type from data set, HDF5 file header attribute 'datatype' may be wrong
00641         sprintf(ipath,"/MDF/images/%d/image",image_index);
00642         hid_t ds=H5Dopen(file,ipath);
00643 
00644         if(ds>0) {      //ds>0 means successfully open the dataset
00645                 hid_t dt = H5Dget_type(ds);
00646 
00647                 switch(H5Tget_size(dt)) {
00648                 case 4:
00649                         dict["datatype"] = (int)EMUtil::EM_FLOAT;
00650                         break;
00651                 case 2:
00652                         dict["datatype"] = (int)EMUtil::EM_USHORT;
00653                         break;
00654                 case 1:
00655                         dict["datatype"] = (int)EMUtil::EM_UCHAR;
00656                         break;
00657                 default:
00658                         throw ImageReadException(filename, "EMAN does not support this data type.");
00659                 }
00660 
00661                 H5Tclose(dt);
00662         }
00663 
00664         H5Dclose(ds);
00665 
00666         EXITFUNC;
00667         return 0;
00668 }
00669 
00670 // This erases any existing attributes from the image group
00671 // prior to writing a new header. For a new image there
00672 // won't be any, so this should be harmless.
00673 int HdfIO2::erase_header(int image_index)
00674 {
00675         ENTERFUNC;
00676 
00677         if(image_index < 0) return 0; //image_index<0 for appending image, no need for erasing
00678 
00679         init();
00680 #ifdef DEBUGHDF
00681         printf("erase_head %d\n",image_index);
00682 #endif
00683         int i;
00684         // Each image is in a group for later expansion. Open the group
00685         char ipath[50];
00686         sprintf(ipath,"/MDF/images/%d", image_index);
00687         hid_t igrp=H5Gopen(file, ipath);
00688 
00689         int nattr=H5Aget_num_attrs(igrp);
00690 
00691         char name[ATTR_NAME_LEN];
00692         for (i=0; i<nattr; i++) {
00693                 hid_t attr = H5Aopen_idx(igrp, 0); //use 0 as index here, since the H5Adelete() will change the index
00694                 H5Aget_name(attr,127,name);
00695                 H5Aclose(attr);
00696                 if( H5Adelete(igrp,name) < 0 ) {
00697                         LOGERR("attribute %s deletion error in erase_header().\n", name);
00698                 }
00699         }
00700 
00701         H5Gclose(igrp);
00702         EXITFUNC;
00703         return 0;
00704 }
00705 
00706 
00707 int HdfIO2::read_data(float *data, int image_index, const Region *area, bool)
00708 {
00709         ENTERFUNC;
00710 #ifdef DEBUGHDF
00711         printf("read_data %d\n",image_index);
00712 #endif
00713 
00714         char ipath[50];
00715         sprintf(ipath,"/MDF/images/%d/image",image_index);
00716         hid_t ds=H5Dopen(file,ipath);
00717         if (ds<0) throw ImageWriteException(filename,"Image does not exist");
00718         hid_t spc=H5Dget_space(ds);
00719         hid_t dt = H5Dget_type(ds);
00720 
00721         hsize_t dims_out[3];
00722         hsize_t rank = H5Sget_simple_extent_ndims(spc);
00723 
00724         H5Sget_simple_extent_dims(spc, dims_out, NULL);
00725         if(rank == 1) {
00726                 nx = dims_out[0];
00727                 ny = 1;
00728                 nz = 1;
00729         }
00730         else if(rank == 2) {
00731                 nx = dims_out[1];
00732                 ny = dims_out[0];
00733                 nz = 1;
00734         }
00735         else if(rank == 3) {
00736                 nx = dims_out[2];
00737                 ny = dims_out[1];
00738                 nz = dims_out[0];
00739         }
00740 
00741         if (area) {
00742                 hid_t memoryspace = 0;
00743 
00744                 /*Get the file dataspace - the region we want to read in the file*/
00745                 int x0 = 0, y0 = 0, z0 = 0;             //the coordinates for up left corner, trim to be within image bound
00746                 int x1 = 0, y1 = 0, z1 = 0;             //the coordinates for down right corner, trim to be within image bound
00747                 int nx1 = 1, ny1 = 1, nz1 = 1;  //dimensions of the sub-region, actual region read from file
00748                 if(rank == 3) {
00749                         hsize_t     doffset[3];             /* hyperslab offset in the file */
00750                         doffset[2] = (hsize_t)(area->x_origin() < 0 ? 0 : area->x_origin());
00751                         doffset[1] = (hsize_t)(area->y_origin() < 0 ? 0 : area->y_origin());
00752                         doffset[0] = (hsize_t)(area->z_origin() < 0 ? 0 : area->z_origin());
00753                         x0 = (int)doffset[0];
00754                         y0 = (int)doffset[1];
00755                         z0 = (int)doffset[2];
00756 
00757                         z1 = (int)(area->x_origin() + area->get_width());
00758                         z1 = (int)(z1 > static_cast<int>(nx) ? nx : z1);
00759 
00760                         y1 = (int)(area->y_origin() + area->get_height());
00761                         y1 = (int)(y1 > static_cast<int>(ny) ? ny : y1);
00762                         if(y1 <= 0) {
00763                                 y1 = 1;
00764                         }
00765 
00766                         x1 = (int)(area->z_origin() + area->get_depth());
00767                         x1 = (int)(x1 > static_cast<int>(nz) ? nz : x1);
00768                         if(x1 <= 0) {
00769                                 x1 = 1;
00770                         }
00771 
00772                         if(x1 < x0 || y1< y0 || z1 < z0) return 0; //out of bounds, this is fine, nothing happens
00773 
00774                         hsize_t     dcount[3];              /* size of the hyperslab in the file */
00775                         dcount[0] = x1 - doffset[0];
00776                         dcount[1] = y1 - doffset[1];
00777                         dcount[2] = z1 - doffset[2];
00778 
00779                         H5Sselect_hyperslab (spc, H5S_SELECT_SET, (const hsize_t*)doffset, NULL, (const hsize_t*)dcount, NULL);
00780 
00781                         /*Define memory dataspace - the memory we will created for the region*/
00782                         hsize_t     dims[3];              /* size of the region in the memory */
00783                         dims[0] = dcount[2]?dcount[2]:1;
00784                         dims[1] = dcount[1]?dcount[1]:1;
00785                         dims[2] = dcount[0]?dcount[0]:1;
00786                         nx1 = (int)dims[0];
00787                         ny1 = (int)dims[1];
00788                         nz1 = (int)dims[2];
00789 
00790                         memoryspace = H5Screate_simple(3, dims, NULL);
00791                 }
00792                 else if(rank == 2) {
00793                         hsize_t     doffset[2];             /* hyperslab offset in the file */
00794                         doffset[1] = (hsize_t)(area->x_origin() < 0 ? 0 : area->x_origin());
00795                         doffset[0] = (hsize_t)(area->y_origin() < 0 ? 0 : area->y_origin());
00796                         x0 = (int)doffset[0];
00797                         y0 = (int)doffset[1];
00798                         z0 = 1;
00799 
00800                         y1 = (int)(area->x_origin() + area->get_width());
00801                         y1 = (int)(y1 > static_cast<int>(nx) ? nx : y1);
00802 
00803                         x1 = (int)(area->y_origin() + area->get_height());
00804                         x1 = (int)(x1 > static_cast<int>(ny) ? ny : x1);
00805                         if(x1 <= 0) {
00806                                 x1 = 1;
00807                         }
00808 
00809                         z1 = 1;
00810 
00811                         if(x1 < x0 || y1< y0) return 0; //out of bounds, this is fine, nothing happens
00812 
00813                         hsize_t     dcount[2];              /* size of the hyperslab in the file */
00814                         dcount[0] = x1 - doffset[0];
00815                         dcount[1] = y1 - doffset[1];
00816 
00817                         H5Sselect_none(spc);
00818                         H5Sselect_hyperslab (spc, H5S_SELECT_SET, (const hsize_t*)doffset, NULL, (const hsize_t*)dcount, NULL);
00819 
00820                         /*Define memory dataspace - the memory we will created for the region*/
00821                         hsize_t     dims[2];              /* size of the region in the memory */
00822                         dims[0] = (hsize_t)(dcount[1]?dcount[1]:1);
00823                         dims[1] = (hsize_t)(dcount[0]?dcount[0]:1);
00824                         nx1 = (int)dims[0];
00825                         ny1 = (int)dims[1];
00826                         nz1 = 1;
00827 
00828                         memoryspace = H5Screate_simple(2, dims, NULL);
00829                 }
00830 
00831                 if( (area->x_origin()>=0) && (area->y_origin()>=0) && (area->z_origin()>=0) && ((hsize_t)(area->x_origin() + area->get_width())<=nx) && ((hsize_t)(area->y_origin() + area->get_height())<=ny) && ((hsize_t)(area->z_origin() + area->get_depth())<=nz) ){      //the region is in boundary
00832                         H5Dread(ds,H5T_NATIVE_FLOAT,memoryspace,spc,H5P_DEFAULT,data);
00833                 }
00834                 else {  //the region are partial out of boundary
00835                         /* When the requested region has some part out of image boundary,
00836                          * we need read the sub-area which is within image,
00837                          * and fill the out of boundary part with zero.
00838                          * We actually read the sub-region from HDF by hyperslab I/O,
00839                          * then copy it back to the pre-allocated region.*/
00840                         float * subdata = new float[nx1*ny1*nz1];
00841 
00842 
00843                         H5Dread(ds,H5T_NATIVE_FLOAT,memoryspace,spc,H5P_DEFAULT,subdata);
00844 
00845                         int xd0=0, yd0=0, zd0=0;        //The coordinates of the top-left corner sub-region in region
00846                         size_t clipped_row_size = 0;
00847                         if(rank == 3) {
00848                                 xd0 = (int) (area->x_origin() < 0 ? -area->x_origin() : 0);
00849                                 yd0 = (int) (area->y_origin() < 0 ? -area->y_origin() : 0);
00850                                 zd0 = (int) (area->z_origin() < 0 ? -area->z_origin() : 0);
00851                                 clipped_row_size = (z1-z0)* sizeof(float);
00852                         }
00853                         else if(rank == 2) {
00854                                 xd0 = (int) (area->x_origin() < 0 ? -area->x_origin() : 0);
00855                                 yd0 = (int) (area->y_origin() < 0 ? -area->y_origin() : 0);
00856                                 clipped_row_size = (y1-y0)* sizeof(float);
00857                         }
00858 
00859                         int src_secsize = nx1 * ny1;
00860                         int dst_secsize = (int)(area->get_width())*(int)(area->get_height());
00861 
00862                         float * src = subdata;
00863                         float * dst = data + zd0*dst_secsize + yd0*(int)(area->get_width()) + xd0;
00864 
00865                         int src_gap = src_secsize - (y1-y0) * nx1;
00866                         int dst_gap = dst_secsize - (y1-y0) * (int)(area->get_width());
00867 
00868                         for(int i = 0; i<nz1; ++i) {
00869                                 for(int j = 0; j<ny1; ++j) {
00870                                         EMUtil::em_memcpy(dst, src, clipped_row_size);
00871 
00872                                         src += nx1;
00873                                         dst += (int)(area->get_width());
00874                                 }
00875                                 src += src_gap;
00876                                 dst += dst_gap;
00877                         }
00878 
00879                         delete [] subdata;
00880                 }
00881                 H5Sclose(memoryspace);
00882         } else {
00883                 hsize_t size = (hsize_t)nx*ny*nz;
00884                 hsize_t i=0;
00885                 hsize_t j=0;
00886                 unsigned short *usdata = (unsigned short *) data;
00887                 unsigned char *cdata = (unsigned char *) data;
00888                 switch(H5Tget_size(dt)) {
00889                 case 4:
00890                         H5Dread(ds,H5T_NATIVE_FLOAT,spc,spc,H5P_DEFAULT,data);
00891                         break;
00892                 case 2:
00893                         H5Dread(ds,H5T_NATIVE_USHORT,spc,spc,H5P_DEFAULT,usdata);
00894                         for (i = 0; i < size; ++i) {
00895                                 j = size - 1 - i;
00896                                 data[j] = static_cast < float >(usdata[j]);
00897                         }
00898                         break;
00899                 case 1:
00900                         H5Dread(ds,H5T_NATIVE_UCHAR,spc,spc,H5P_DEFAULT,cdata);
00901                         for (i = 0; i < size; ++i) {
00902                                 j = size - 1 - i;
00903                                 data[j] = static_cast < float >(cdata[j]);
00904                         }
00905                         break;
00906                 default:
00907                         throw ImageReadException(filename, "EMAN does not support this data type.");
00908                 }
00909         }
00910 
00911         H5Tclose(dt);
00912         H5Sclose(spc);
00913         H5Dclose(ds);
00914         EXITFUNC;
00915         return 0;
00916 }
00917 
00918 
00919 // Writes all attributes in 'dict' to the image group
00920 // Creation of the image dataset is also handled here
00921 int HdfIO2::write_header(const Dict & dict, int image_index, const Region* area,
00922                                                 EMUtil::EMDataType, bool)
00923 {
00924 #ifdef DEBUGHDF
00925         printf("write_head %d\n",image_index);
00926 #endif
00927         ENTERFUNC;
00928         init();
00929 
00930         nx = (int)dict["nx"];
00931         ny = (int)dict["ny"];
00932         nz = (int)dict["nz"];
00933 
00934         if(image_index<0) {
00935                 image_index = get_nimg();
00936         }
00937 
00938         // If image_index<0 append, and make sure the max value in the file is correct
00939         // though this is normally handled by EMData.write_image()
00940         hid_t attr=H5Aopen_name(group,"imageid_max");
00941         int nimg = read_attr(attr);
00942         H5Aclose(attr);
00943 
00944         unsigned int i;
00945         if (image_index<0) image_index=nimg+1;
00946         if (image_index>nimg) {
00947                 write_attr(group,(const char *)"imageid_max",EMObject(image_index));
00948         }
00949 
00950         // Each image is in a group for later expansion. Open the group, create if necessary
00951         char ipath[50];
00952         sprintf(ipath,"/MDF/images/%d",image_index);
00953         hid_t igrp=H5Gopen(file,ipath);
00954 
00955         if (igrp<0) {   //group not existed
00956                 is_exist = false;
00957                 // Need to create a new image group
00958                 igrp=H5Gcreate(file,ipath,64);          // The image is a group, with attributes on the group
00959                 if (igrp<0) throw ImageWriteException(filename,"Unable to add /MDF/images/# to HDF5 file");
00960         }
00963         else {
00964                 is_exist = true;
00965                 int nattr=H5Aget_num_attrs(igrp);
00966                 char name[ATTR_NAME_LEN];
00967                 Dict dict2;
00968                 for (int i=0; i<nattr; i++) {
00969                         hid_t attr=H5Aopen_idx(igrp, i);
00970                         H5Aget_name(attr,127,name);
00971                         if (strncmp(name,"EMAN.", 5)!=0) {
00972                                 H5Aclose(attr);
00973                                 continue;
00974                         }
00975                         EMObject val=read_attr(attr);
00976                         dict2[name+5]=val;
00977                         H5Aclose(attr);
00978 
00979                         if(!dict2.has_key("datatype")) {//by default, HDF5 is written as float
00980                                 dict2["datatype"] = (int)EMUtil::EM_FLOAT;
00981                         }
00982                 }
00983 
00984                 if(area) {
00985                         check_region(area, IntSize(dict2["nx"], dict2["ny"], dict2["nz"]), false, true);
00986                 }
00987                 else {
00988                         erase_header(image_index);
00989 
00990                         //change the size or data type of a image,
00991                         //the existing data set is invalid, unlink it
00992                         if( (int)dict["nx"]*(int)dict["ny"]*(int)dict["nz"] !=
00993                                 (int)dict2["nx"]*(int)dict2["ny"]*(int)dict2["nz"] ||
00994                                 dict["datatype"] != dict2["datatype"] ) {
00995                                 sprintf(ipath,"/MDF/images/%d/image",image_index);
00996                                 H5Gunlink(igrp, ipath);
00997                         }
00998                 }
00999         }
01000 
01001         if(!area) {
01002                 // Write the attributes to the group
01003                 vector <string> keys=dict.keys();
01004 
01005                 for (i=0; i<keys.size(); i++) {
01006                         string s("EMAN.");
01007                         s+=keys[i];
01008                         write_attr(igrp,s.c_str(),dict[keys[i]]);
01009                 }
01010         }
01011 
01012         H5Gclose(igrp);
01013 
01014     // Set render_min and render_max from EMData attr's if possible.
01015         if(dict.has_key("render_min")) rendermin=(float)dict["render_min"];
01016         else rendermin=0;
01017         if(dict.has_key("render_max")) rendermax=(float)dict["render_max"];
01018         else rendermax=0;
01019 
01020         EXITFUNC;
01021         return 0;
01022 }
01023 
01024 // Writes the actual image data to the corresponding dataset (already created)
01025 int HdfIO2::write_data(float *data, int image_index, const Region* area,
01026                                           EMUtil::EMDataType dt, bool)
01027 {
01028         ENTERFUNC;
01029         if(!data) {
01030                 std::cerr << "Warning:blank image written!!! " << std::endl;
01031                 return 0;
01032         }
01033 
01034 #ifdef DEBUGHDF
01035         printf("write_data %d\n",image_index);
01036 #endif
01037 
01038         if (image_index<0) {
01039                 hid_t attr=H5Aopen_name(group,"imageid_max");
01040                 image_index = read_attr(attr);
01041                 H5Aclose(attr);
01042         }
01043 
01044         hid_t spc;      //dataspace
01045         hid_t ds;       //dataset
01046         char ipath[50];
01047         sprintf(ipath,"/MDF/images/%d/image",image_index);
01048 
01049         // Now create the actual image dataspace(not for regional writing)
01050         if (nz==1)  {
01051                 hsize_t dims[2]= { ny,nx };
01052                 spc=H5Screate_simple(2,dims,NULL);
01053         }
01054         else {
01055                 hsize_t dims[3]= { nz, ny, nx };
01056                 spc=H5Screate_simple(3,dims,NULL);
01057         }
01058 
01059         ds=H5Dopen(file,ipath);
01060         hsize_t rank = 0;
01061         if(ds<0) {//new dataset
01062                 switch(dt) {
01063                 case EMUtil::EM_FLOAT:
01064                         ds=H5Dcreate(file,ipath, H5T_NATIVE_FLOAT, spc, H5P_DEFAULT );
01065                         break;
01066                 case EMUtil::EM_USHORT:
01067                         ds=H5Dcreate(file,ipath, H5T_NATIVE_USHORT, spc, H5P_DEFAULT );
01068                         break;
01069                 case EMUtil::EM_UCHAR:
01070                         ds=H5Dcreate(file,ipath, H5T_NATIVE_UCHAR, spc, H5P_DEFAULT );
01071                         break;
01072                 default:
01073                         throw ImageWriteException(filename,"HDF5 does not support this data format");
01074                 }
01075         }
01076         else {  //existing file
01077                 hid_t spc_file = H5Dget_space(ds);
01078                 rank = H5Sget_simple_extent_ndims(spc_file);
01079                 H5Sclose(spc_file);
01080         }
01081 
01082         //convert data to unsigned short, unsigned char...
01083         hsize_t size = (hsize_t)nx*ny*nz;
01084         unsigned char *cdata = 0;
01085         unsigned short *usdata = 0;
01086         // float rendermin = 0.0f;
01087         // float rendermax = 0.0f;
01088         if (!rendermin && !rendermax) EMUtil::getRenderMinMax(data, nx, ny, rendermin, rendermax, nz);
01089 
01090         if(area) {
01091                 hid_t filespace = H5Dget_space(ds);
01092                 hid_t memoryspace = 0;
01093                 if(rank==3) {
01094                         hsize_t doffset[3];             /*hyperslab offset in the file*/
01095                         doffset[0] = (hsize_t)(area->z_origin()<0 ? 0 : area->z_origin());
01096                         doffset[1] = (hsize_t)(area->y_origin()<0 ? 0 : area->y_origin());
01097                         doffset[2] = (hsize_t)(area->x_origin()<0 ? 0 : area->x_origin());
01098 
01099                         hsize_t dcount[3];              /*size of the hyperslab in the file*/
01100                         dcount[0] = (hsize_t)(area->get_depth()?area->get_depth():1);
01101                         dcount[1] = (hsize_t)(area->get_height()?area->get_height():1);
01102                         dcount[2] = (hsize_t)(area->get_width()?area->get_width():1);
01103 
01104                         herr_t err_no = H5Sselect_hyperslab(filespace, H5S_SELECT_SET, (const hsize_t*)doffset, NULL, (const hsize_t*)dcount, NULL);
01105                         if(err_no < 0) {
01106                                 std::cerr << "H5Sselect_hyperslab error: " << err_no << std::endl;
01107                         }
01108 
01109                         /*Create memory space with size of the region.*/
01110                         hsize_t dims[3];        /*size of the region in the memory*/
01111                         dims[0] = dcount[2]?dcount[2]:1;
01112                         dims[1] = dcount[1]?dcount[1]:1;
01113                         dims[2] = dcount[0]?dcount[0]:1;
01114 
01115                         memoryspace = H5Screate_simple(rank, dims, NULL);
01116                 }
01117                 else if(rank==2){
01118                         hsize_t doffset[2];             /*hyperslab offset in the file*/
01119                         doffset[0] = (hsize_t)(area->y_origin() < 0 ? 0 : area->y_origin());
01120                         doffset[1] = (hsize_t)(area->x_origin() < 0 ? 0 : area->x_origin());
01121 
01122                         hsize_t dcount[2];              /*size of the hyperslab in the file*/
01123                         dcount[0] = (hsize_t)area->get_height();
01124                         dcount[1] = (hsize_t)area->get_width();
01125 
01126                         herr_t err_no = H5Sselect_hyperslab(filespace, H5S_SELECT_SET, (const hsize_t*)doffset, NULL, (const hsize_t*)dcount, NULL);
01127                         if(err_no < 0) {
01128                                 std::cerr << "H5Sselect_hyperslab error: " << err_no << std::endl;
01129                         }
01130 
01131                         /*Create memory space with size of the region.*/
01132                         /*Define memory dataspace - the memory we will created for the region*/
01133                         hsize_t     dims[2];              /* size of the region in the memory */
01134                         dims[0] = (hsize_t)(dcount[1]?dcount[1]:1);
01135                         dims[1] = (hsize_t)(dcount[0]?dcount[0]:1);
01136 
01137                         memoryspace = H5Screate_simple(rank, dims, NULL);
01138                 }
01139                 else {
01140                         std::cerr << "rank is wrong: " << rank << std::endl;
01141                 }
01142 
01143                 herr_t err_no;
01144                 switch(dt) {
01145                 case EMUtil::EM_FLOAT:
01146                         err_no = H5Dwrite(ds, H5T_NATIVE_FLOAT, memoryspace, filespace, H5P_DEFAULT, data);
01147                         if(err_no<0) {
01148                                 std::cerr << "H5Dwrite error: " << err_no << std::endl;
01149                         }
01150                         break;
01151                 default:
01152                         throw ImageWriteException(filename,"HDF5 does not support regional writing for this data format");
01153                 }
01154                 H5Sclose(filespace);
01155                 H5Sclose(memoryspace);
01156         }
01157         else {
01158                 switch(dt) {
01159                 case EMUtil::EM_FLOAT:
01160                         H5Dwrite(ds,H5T_NATIVE_FLOAT,spc,spc,H5P_DEFAULT,data);
01161                         break;
01162                 case EMUtil::EM_USHORT:
01163                         usdata = new unsigned short[size];
01164                         for (size_t i = 0; i < size; ++i) {
01165                                 if(data[i] <= rendermin) {
01166                                         usdata[i] = 0;
01167                                 }
01168                                 else if(data[i] >= rendermax) {
01169                                         usdata[i] = USHRT_MAX;
01170                                 }
01171                                 else {
01172                                         usdata[i]=(unsigned short)((data[i]-rendermin)/(rendermax-rendermin)*USHRT_MAX);
01173                                 }
01174                         }
01175                         H5Dwrite(ds,H5T_NATIVE_USHORT,spc,spc,H5P_DEFAULT,usdata);
01176                         if(usdata) {delete [] usdata; usdata=0;}
01177                         break;
01178                 case EMUtil::EM_UCHAR:
01179                         cdata = new unsigned char[size];
01180                         for (size_t i = 0; i < size; ++i) {
01181                                 if(data[i] <= rendermin) {
01182                                         cdata[i] = 0;
01183                                 }
01184                                 else if(data[i] >= rendermax){
01185                                         cdata[i] = UCHAR_MAX;
01186                                 }
01187                                 else {
01188                                         cdata[i]=(unsigned char)((data[i]-rendermin)/(rendermax-rendermin)*UCHAR_MAX);
01189                                 }
01190                         }
01191                         H5Dwrite(ds,H5T_NATIVE_UCHAR,spc,spc,H5P_DEFAULT,cdata);
01192                         if(cdata) {delete [] cdata; cdata=0;}
01193                         break;
01194                 default:
01195                         throw ImageWriteException(filename,"HDF5 does not support this data format");
01196                 }
01197         }
01198 
01199         H5Sclose(spc);
01200         H5Dclose(ds);
01201         EXITFUNC;
01202         return 0;
01203 }
01204 
01205 int HdfIO2::get_nimg()
01206 {
01207         init();
01208         hid_t attr=H5Aopen_name(group,"imageid_max");
01209         int n = read_attr(attr);
01210         H5Aclose(attr);
01211 
01212         return n+1;
01213 }
01214 
01215 void HdfIO2::flush()
01216 {
01217         return;
01218 }
01219 
01220 bool HdfIO2::is_complex_mode()
01221 {
01222         return false;
01223 }
01224 
01225 // always big endian
01226 bool HdfIO2::is_image_big_endian()
01227 {
01228         return true;
01229 }
01230 
01231 
01232 
01233 #endif  //EM_HDF5

Generated on Tue Jun 11 13:40:38 2013 for EMAN2 by  doxygen 1.3.9.1