EMAN2
amiraio.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 #include "amiraio.h"
00037 #include "util.h"
00038 
00039 #ifndef WIN32
00040         #include <sys/param.h>
00041 #else
00042         #include <windows.h>
00043         #define MAXPATHLEN (MAX_PATH*4)
00044 #endif  //WIN32
00045 
00046 #include <cstdio>
00047 
00048 using namespace EMAN;
00049 
00050 const char *AmiraIO::MAGIC = "# AmiraMesh";
00051 
00052 AmiraIO::AmiraIO(const string & file, IOMode rw)
00053 :       filename(file), rw_mode(rw), amira_file(0),
00054         is_big_endian(true), initialized(false), dt(EMUtil::EM_UNKNOWN),
00055          nx(0), ny(0), nz(0),
00056          pixel(0), xorigin(0), yorigin(0), zorigin(0)
00057 {
00058 }
00059 
00060 AmiraIO::~AmiraIO()
00061 {
00062         if (amira_file) {
00063                 fclose(amira_file);
00064                 amira_file = 0;
00065         }
00066 }
00067 
00068 void AmiraIO::init()
00069 {
00070         ENTERFUNC;
00071 
00072         if (initialized) {
00073                 return;
00074         }
00075 
00076         initialized = true;
00077         bool is_new_file = false;
00078 
00079         amira_file = sfopen(filename, rw_mode, &is_new_file, true);
00080 
00081         if (!is_new_file) {
00082                 char buf[MAXPATHLEN];
00083                 if (!fgets(buf, MAXPATHLEN, amira_file)) {
00084                         throw ImageReadException(filename, "Amira Header");
00085                 }
00086 
00087                 if (!is_valid(buf)) {
00088                         throw ImageReadException(filename, "invalid Amira Mesh file");
00089                 }
00090 
00091                 if(strstr(buf,"BINARY-LITTLE-ENDIAN")!=0) {
00092                         is_big_endian = false;
00093                 }
00094                 else if(strstr(buf,"BINARY")!=0) {
00095                         is_big_endian = true;
00096                 }
00097                 else if(strstr(buf,"2.0")!=0) {
00098                         is_big_endian = true; //AMIRA convention(<=2.0): BIGENDIAN only
00099                 }
00100         }
00101         EXITFUNC;
00102 }
00103 
00104 bool AmiraIO::is_valid(const void *first_block)
00105 {
00106         ENTERFUNC;
00107         bool result = false;
00108         if (!first_block) {
00109                 result = false;
00110         }
00111         else {
00112                 result = Util::check_file_by_magic(first_block, MAGIC);
00113         }
00114         EXITFUNC;
00115         return result;
00116 }
00117 
00118 int AmiraIO::read_header(Dict & dict, int image_index, const Region *, bool)
00119 {
00120         ENTERFUNC;
00121 
00122         //single image format, index can only be zero
00123         if(image_index == -1) {
00124                 image_index = 0;
00125         }
00126 
00127         if(image_index != 0) {
00128                 throw ImageReadException(filename, "no stack allowed for MRC image. For take 2D slice out of 3D image, read the 3D image first, then use get_clip().");
00129         }
00130 
00131         init();
00132 
00133         char ll[MAXPATHLEN+10];
00134         char datatype[16] = "";
00135 
00136         do {
00137                 fgets(ll,MAXPATHLEN,amira_file);
00138 //              printf("%s", ll);
00139                 if(char* s=strstr(ll,"define Lattice ")) {
00140                         if(sscanf(s+15,"%d %d %d",&nx, &ny, &nz) == 3) {
00141                                 dict["nx"] = nx;
00142                                 dict["ny"] = ny;
00143                                 dict["nz"] = nz;
00144 //                              printf("nx = %d, ny = %d, nz = %d\n", nx, ny, nz);
00145                         };
00146                 }
00147                 else if(char* s=strstr(ll,"BoundingBoxXY")) {           //amiramesh 2.1
00148                         float bx0, bx1, by0, by1;
00149                         if (sscanf(s+13,"%f %f %f %f",&bx0, &bx1, &by0, &by1) == 4 ) {
00150                                 pixel = (bx1-bx0)/(nx-1);
00151                                 xorigin = bx0;
00152                                 yorigin = by0;
00153                                 dict["apix_x"] = pixel;
00154                                 dict["apix_y"] = pixel;
00155                                 dict["apix_z"] = pixel;
00156                                 dict["origin_x"] = xorigin;
00157                                 dict["origin_y"] = yorigin;
00158 //                              printf("pixel=%g: xorigin=%g yorigin=%g \n", pixel, xorigin, yorigin);
00159                         }
00160                 }
00161                 else if(char* s=strstr(ll,"BoundingBox")) {             //amiramesh 2.0
00162                         float bx0, bx1, by0, by1, bz0, bz1;
00163                         if (sscanf(s+11,"%f %f %f %f %f %f",&bx0, &bx1, &by0, &by1, &bz0, &bz1) == 6 ) {
00164                                 pixel = (bx1-bx0)/(nx-1);
00165                                 xorigin = bx0;
00166                                 yorigin = by0;
00167                                 zorigin = bz0;
00168                                 dict["apix_x"] = pixel;
00169                                 dict["apix_y"] = pixel;
00170                                 dict["apix_z"] = pixel;
00171                                 dict["origin_x"] = xorigin;
00172                                 dict["origin_y"] = yorigin;
00173                                 dict["origin_z"] = zorigin;
00174 //                              printf("pixel=%g: xorigin=%g yorigin=%g zorigin=%g\n", pixel, xorigin, yorigin, zorigin);
00175                         }
00176                 }
00177                 else if(char* s=strstr(ll,"Lattice { ")) {
00178                         sscanf(s+10,"%s ", datatype);
00179                         if(!strncmp(datatype, "float", 5)) {
00180                                 dt = EMUtil::EM_FLOAT;
00181                                 dict["datatype"] = EMUtil::EM_FLOAT;
00182                         }
00183                         else if(!strncmp(datatype, "short", 5)) {
00184                                 dt = EMUtil::EM_SHORT;
00185                                 dict["datatype"] = EMUtil::EM_SHORT;
00186                         }
00187                         else if(!strncmp(datatype, "byte", 4)) {
00188                                 dt = EMUtil::EM_CHAR;
00189                                 dict["datatype"] = EMUtil::EM_CHAR;
00190                         }
00191                         else {
00192                                 fprintf(stderr,"AmiraIO::read_header: data type \"%s\" is not supported yet\n", datatype);
00193                                 return -1;
00194                         }
00195                 }
00196         } while (! ( ll[0]=='@' && ll[1]=='1') );
00197 
00198         EXITFUNC;
00199         return 0;
00200 
00201 }
00202 
00203 int AmiraIO::write_header(const Dict & dict, int image_index, const Region*, EMUtil::EMDataType, bool)
00204 {
00205         ENTERFUNC;
00206         int err = 0;
00207 
00208         //single image format, index can only be zero
00209         if(image_index == -1) {
00210                 image_index = 0;
00211         }
00212 
00213         if(image_index != 0) {
00214                 throw ImageWriteException(filename, "Amira file does not support stack.");
00215         }
00216         check_write_access(rw_mode, image_index, 1);
00217 
00218         nx = dict["nx"];
00219         ny = dict["ny"];
00220         nz = dict["nz"];
00221 
00222         float xorigin = 0.0f;
00223         if(dict.has_key("origin_x")) xorigin = dict["origin_x"];
00224         float yorigin = 0.0f;
00225         if(dict.has_key("origin_y")) yorigin = dict["origin_y"];
00226         float zorigin = 0.0f;
00227         if(dict.has_key("origin_z")) zorigin = dict["origin_z"];
00228         float pixel = 0.0f;
00229         if(dict.has_key("apix_x")) pixel = dict["apix_x"];
00230 
00231         rewind(amira_file);
00232 
00233         string line1;
00234         if(ByteOrder::is_host_big_endian()) {
00235                 line1= "# AmiraMesh 3D BINARY 2.1\n\n";
00236         }
00237         else {
00238                 line1= "# AmiraMesh BINARY-LITTLE-ENDIAN 2.1\n\n";
00239         }
00240 
00241         string type = "float";
00242 
00243         if (fprintf(amira_file, "%s", line1.c_str()) <= 0) {
00244                 LOGERR("cannot write to AmiraMesh file '%s'", filename.c_str());
00245                 err = 1;
00246         }
00247         else {
00248                 fprintf(amira_file,"define Lattice %d %d %d\n\n",nx,ny,nz);
00249                 fprintf(amira_file,"Parameters {\n");
00250                 fprintf(amira_file, "\tContent \"%dx%dx%d %s, uniform coordinates\",\n", nx,ny,nz,type.c_str());
00251                 fprintf(amira_file, "\tCoordType \"uniform\",\n");
00252                 fprintf(amira_file,"\tBoundingBox %.2f %.2f %.2f %.2f %.2f %.2f\n}\n\n",xorigin,xorigin+pixel*(nx-1),yorigin,yorigin+pixel*(ny-1),zorigin,zorigin+pixel*(nz-1));
00253 
00254                 fprintf(amira_file, "Lattice { float ScalarField } @1\n\n# Data section follows\n@1\n");
00255         }
00256 
00257         EXITFUNC;
00258         return err;
00259 }
00260 
00261 int AmiraIO::read_data(float * rdata, int, const Region *, bool)
00262 {
00263         ENTERFUNC;
00264 
00265         size_t size = (size_t)nx*ny*nz;
00266         switch(dt) {
00267         case EMUtil::EM_FLOAT:
00268                 
00269                 fread(rdata,nx*nz,ny*sizeof(float),amira_file);
00270                 
00271                 if( (is_big_endian && ByteOrder::is_host_big_endian()) && !(is_big_endian || ByteOrder::is_host_big_endian()) ) {
00272                         char tmpdata;
00273                         char *cdata=(char*)rdata;
00274                         for(size_t i=0;i<size;++i){
00275                                 //swap 0-3
00276                                 tmpdata=cdata[4*i+3];
00277                                 cdata[4*i+3]=cdata[4*i];
00278                                 cdata[4*i] = tmpdata;
00279                                 //swap 1-2
00280                                 tmpdata=cdata[4*i+2];
00281                                 cdata[4*i+2]=cdata[4*i+1];
00282                                 cdata[4*i+1] = tmpdata;
00283                                 
00284                                 
00285                         }
00286                 }
00287                 
00288                 float* copydata;
00289                 copydata=new float[size];
00290                 
00291                 for (size_t i=0;i<size;i++){
00292                         copydata[i]=rdata[i];   
00293                 }
00294                 
00295                 for (int i=0;i<ny;i++){
00296                         for (int j=0;j<nx;j++){
00297                                 rdata[i*nx+j]=copydata[(ny-1-i)*nx+j];
00298                                 }
00299                 }
00300                 delete[] copydata;
00301                 
00302                 break;
00303         case EMUtil::EM_SHORT:
00304         {
00305                 short *datashort = (short*)malloc(sizeof(short)*nx*ny*nz);
00306                 fread(datashort,nx*nz,ny*sizeof(short),amira_file);
00307                 if( (is_big_endian && ByteOrder::is_host_big_endian()) && !(is_big_endian || ByteOrder::is_host_big_endian()) ) {
00308                         char tmpdata;
00309                         char *cdata=(char*)datashort;
00310                         for(size_t i=0;i<size;++i){
00311                                 //swap 0-1
00312                                 tmpdata=cdata[2*i+1];
00313                                 cdata[2*i+1]=cdata[2*i];
00314                                 cdata[2*i] = tmpdata;
00315                         }
00316                         for(size_t i=0; i<size; ++i) rdata[i] = float(datashort[i]);
00317                         free(datashort);
00318                 }
00319                 
00320                 
00321                 float* copydata;
00322                 copydata=new float[size];
00323                 
00324                 for (size_t i=0;i<size;i++){
00325                         copydata[i]=rdata[i];   
00326                 }
00327                 
00328                 for (int i=0;i<ny;i++){
00329                         for (int j=0;j<nx;j++){
00330                                 rdata[i*nx+j]=copydata[(ny-1-i)*nx+j];
00331                                 }
00332                 }
00333                 delete[] copydata;
00334                 
00335         }
00336                 break;
00337         case EMUtil::EM_CHAR:
00338         {
00339                 char *databyte = (char*)malloc(sizeof(char)*nx*ny*nz);
00340                 fread(databyte,nx*nz,ny*sizeof(char),amira_file);
00341                 for(size_t i=0; i<size; ++i) rdata[i] = float(databyte[i]);
00342                 free(databyte);
00343                 
00344                 float* copydata;
00345                 copydata=new float[size];
00346                 
00347                 for (size_t i=0;i<size;i++){
00348                         copydata[i]=rdata[i];   
00349                 }
00350                 
00351                 for (int i=0;i<ny;i++){
00352                         for (int j=0;j<nx;j++){
00353                                 rdata[i*nx+j]=copydata[(ny-1-i)*nx+j];
00354                                 }
00355                 }
00356                 delete[] copydata;
00357                 
00358         }
00359                 break;
00360         default:
00361                 fprintf(stderr,"AmiraIO::read_data: data type is not supported yet\n");
00362                 return -1;
00363         }
00364 
00365         EXITFUNC;
00366         return 0;
00367 }
00368 
00369 int AmiraIO::write_data(float *data, int image_index, const Region*, EMUtil::EMDataType, bool)
00370 {
00371         ENTERFUNC;
00372 
00373         check_write_access(rw_mode, image_index, 1, data);
00374 //      ByteOrder::become_big_endian(data, nx * ny * nz);
00375         size_t size = (size_t)nx*ny*nz;
00376         
00377         float* copydata;
00378         copydata=new float[size];
00379         
00380         for (size_t i=0;i<size;i++){
00381                 copydata[i]=data[i];    
00382         }
00383         
00384         for (int i=0;i<ny;i++){
00385                 for (int j=0;j<nx;j++){
00386                         data[i*nx+j]=copydata[(ny-1-i)*nx+j];
00387                         }
00388         }
00389         delete[] copydata;      
00390         
00391         if (fwrite(data, nx * nz, ny * sizeof(float), amira_file) != ny * sizeof(float)) {
00392                 throw ImageWriteException(filename, "incomplete file write in AmiraMesh file");
00393         }
00394 
00395         EXITFUNC;
00396         return 0;
00397 }
00398 
00399 void AmiraIO::flush()
00400 {
00401         fflush(amira_file);
00402 }
00403 
00404 bool AmiraIO::is_complex_mode()
00405 {
00406         return false;
00407 }
00408 
00409 bool AmiraIO::is_image_big_endian()
00410 {
00411         init();
00412         return is_big_endian;
00413 }
00414