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

pifio.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 "pifio.h"
00037 #include "portable_fileio.h"
00038 #include "geometry.h"
00039 #include "imageio.h"
00040 #include <cstring>
00041 
00042 #ifdef WIN32
00043 #include <time.h>
00044 #endif
00045 
00046 using namespace EMAN;
00047 
00048 PifIO::PifIO(const string & pif_filename, IOMode rw)
00049 :       filename(pif_filename), rw_mode(rw)
00050 {
00051         pif_file = 0;
00052         mode_size = 0;
00053         is_big_endian = ByteOrder::is_host_big_endian();
00054         initialized = false;
00055         real_scale_factor = 1;
00056         is_new_file = false;
00057         memset(&pfh, 0, sizeof(PifFileHeader));
00058 }
00059 
00060 PifIO::~PifIO()
00061 {
00062         if (pif_file) {
00063                 fclose(pif_file);
00064                 pif_file = 0;
00065         }
00066 }
00067 
00068 
00069 int PifIO::get_mode_size(PifDataMode mode)
00070 {
00071         int size = 0;
00072 
00073         switch (mode) {
00074         case PIF_CHAR:
00075         case PIF_BOXED_DATA:
00076                 size = sizeof(char);
00077                 break;
00078         case PIF_SHORT:
00079         case PIF_SHORT_FLOAT:
00080         case PIF_SHORT_COMPLEX:
00081         case PIF_SHORT_FLOAT_COMPLEX:
00082         case PIF_MAP_FLOAT_SHORT:
00083                 size = sizeof(short);
00084                 break;
00085         case PIF_FLOAT:
00086         case PIF_FLOAT_COMPLEX:
00087                 size = sizeof(float);
00088                 break;
00089         case PIF_FLOAT_INT:
00090         case PIF_FLOAT_INT_COMPLEX:
00091         case PIF_MAP_FLOAT_INT:
00092                 size = sizeof(int);
00093                 break;
00094         default:
00095                 break;
00096         }
00097         return size;
00098 }
00099 
00100 bool PifIO::is_float_int(int m)
00101 {
00102         PifDataMode mode = static_cast < PifDataMode > (m);
00103         switch (mode) {
00104         case PIF_SHORT_FLOAT:
00105         case PIF_SHORT_FLOAT_COMPLEX:
00106         //case PIF_FLOAT:
00107         case PIF_FLOAT_INT:
00108         //case PIF_FLOAT_COMPLEX:
00109         case PIF_FLOAT_INT_COMPLEX:
00110         case PIF_MAP_FLOAT_SHORT:
00111         case PIF_MAP_FLOAT_INT:
00112                 return true;
00113         default:
00114                 break;
00115         }
00116         return false;
00117 }
00118 
00119 void PifIO::init()
00120 {
00121         ENTERFUNC;
00122         if (initialized) {
00123                 return;
00124         }
00125 
00126         initialized = true;
00127         pif_file = sfopen(filename, rw_mode, &is_new_file);
00128 
00129         if (!is_new_file) {
00130                 if (fread(&pfh, sizeof(PifFileHeader), 1, pif_file) != 1) {
00131                         throw ImageReadException(filename, "PIF file header");
00132                 }
00133 
00134                 if (!is_valid(&pfh)) {
00135                         throw ImageReadException(filename, "invalid PIF file");
00136                 }
00137 
00138                 is_big_endian = ByteOrder::is_data_big_endian(&pfh.nz);
00139                 become_host_endian(&pfh.htype);
00140 
00141                 if (pfh.htype != 1) {
00142                         string desc = "only support PIF with all projects having the same dimensions";
00143                         throw ImageReadException(filename, desc);
00144                 }
00145 
00146                 become_host_endian(&pfh.mode);
00147                 become_host_endian(&pfh.nx);
00148                 become_host_endian(&pfh.ny);
00149                 become_host_endian(&pfh.nz);
00150                 become_host_endian(&pfh.nimg);
00151 
00152                 if (is_float_int(pfh.mode)) {
00153                         real_scale_factor = (float) atof(pfh.scalefactor);
00154                 }
00155 
00156                 mode_size = get_mode_size(static_cast < PifDataMode > (pfh.mode));
00157 
00158                 if (is_complex_mode()) {
00159                         pfh.nx *= 2;
00160                 }
00161         }
00162         EXITFUNC;
00163 }
00164 
00165 bool PifIO::is_valid(const void *first_block)
00166 {
00167         ENTERFUNC;
00168         bool result = false;
00169 
00170         if (first_block) {
00171                 const int *data = static_cast < const int *>(first_block);
00172                 int m1 = data[0];
00173                 int m2 = data[1];
00174                 int endian = data[7];
00175                 bool data_big_endian = false;
00176                 if (endian) {
00177                         data_big_endian = true;
00178                 }
00179 
00180                 if (data_big_endian != ByteOrder::is_host_big_endian()) {
00181                         ByteOrder::swap_bytes(&m1);
00182                         ByteOrder::swap_bytes(&m2);
00183                 }
00184 
00185                 if (m1 == PIF_MAGIC_NUM && m2 == PIF_MAGIC_NUM) {
00186                          result = true;
00187                 }
00188         }
00189 
00190         EXITFUNC;
00191         return result;
00192 }
00193 
00194 void PifIO::fseek_to(int image_index)
00195 {
00196         int pih_sz = sizeof(PifImageHeader);
00197         int image_size = 0;
00198 
00199 #if 0
00200         // this works for some images that PURDUE people gave to me.
00201         // But those images don't follow the PIF specification. So
00202         // I believe they are in wrong format.
00203         if (pfh.nimg == 1) {
00204                 image_size = pfh.nx * pfh.ny * pfh.nz;
00205         }
00206         else {
00207                 image_size = pfh.nx * pfh.ny;
00208         }
00209 #endif
00210         image_size = pfh.nx * pfh.ny * pfh.nz;
00211 
00212         size_t file_offset = sizeof(PifFileHeader) +
00213                 (pih_sz + image_size * mode_size) * image_index;
00214 
00215         portable_fseek(pif_file, file_offset, SEEK_SET);
00216 }
00217 
00218 
00219 int PifIO::read_header(Dict & dict, int image_index, const Region * area, bool)
00220 {
00221         ENTERFUNC;
00222 
00223         check_read_access(image_index);
00224         fseek_to(image_index);
00225 
00226         int pih_sz = sizeof(PifImageHeader);
00227         PifImageHeader pih;
00228 
00229         if (fread(&pih, pih_sz, 1, pif_file) != 1) {
00230                 throw ImageReadException(filename, "PIF Image header");
00231         }
00232         else {
00233                 check_region(area, FloatSize(pih.nx, pih.ny, pih.nz), is_new_file);
00234                 int xlen = 0, ylen = 0, zlen = 0;
00235                 EMUtil::get_region_dims(area, pih.nx, &xlen, pih.ny, &ylen, pih.nz, &zlen);
00236 
00237                 dict["nx"] = xlen;
00238                 dict["ny"] = ylen;
00239                 dict["nz"] = zlen;
00240 
00241                 dict["datatype"] = to_em_datatype(pih.mode);
00242 
00243                 dict["apix_x"] = static_cast < float >(pih.xlen);
00244                 dict["apix_y"] = static_cast < float >(pih.ylen);
00245                 dict["apix_z"] = static_cast < float >(pih.zlen);
00246 
00247                 dict["minimum"] = static_cast < float >(pih.min);
00248                 dict["maximum"] = static_cast < float >(pih.max);
00249                 dict["mean"] = static_cast < float >(pih.mean);
00250                 dict["sigma"] = static_cast < float >(pih.sigma);
00251 
00252                 dict["origin_x"] = static_cast < float >(pih.xorigin);
00253                 dict["origin_y"] = static_cast < float >(pih.yorigin);
00254         }
00255 
00256         EXITFUNC;
00257 
00258         return 0;
00259 }
00260 
00261 int PifIO::write_header(const Dict & dict, int image_index, const Region* area,
00262                                                 EMUtil::EMDataType, bool)
00263 {
00264         ENTERFUNC;
00265 
00266         if(image_index<0) {
00267                 image_index = get_nimg();
00268         }
00269 
00270         check_write_access(rw_mode, image_index);
00271 
00272         if (area) {
00273                 check_region(area, FloatSize(pfh.nx, pfh.ny, pfh.nz), is_new_file);
00274                 EXITFUNC;
00275                 return 0;
00276         }
00277 
00278         time_t t0 = time(0);
00279         struct tm *t = localtime(&t0);
00280 
00281         if (!is_new_file) {
00282                 if (is_big_endian != ByteOrder::is_host_big_endian()) {
00283                         throw ImageWriteException(filename, "writing to opposite byteorder file");
00284                 }
00285 
00286                 int nx1 = dict["nx"];
00287                 int ny1 = dict["ny"];
00288                 int nz1 = dict["nz"];
00289 
00290                 int mode1 = to_pif_datatype(dict["datatype"]);
00291 
00292                 if (nx1 != pfh.nx || ny1 != pfh.ny || nz1 != pfh.nz) {
00293                         LOGERR("cannot write to different size file %s (%dx%dx%d != %dx%dx%d)",
00294                                    filename.c_str(), nx1, ny1, nz1, pfh.nx, pfh.ny, pfh.nz);
00295                         throw ImageWriteException(filename, "different image size");
00296                 }
00297 
00298                 if (mode1 != pfh.mode) {
00299                         throw ImageWriteException(filename, "different data type");
00300                 }
00301 
00302                 fseek_to(image_index);
00303         }
00304         else {
00305                 pfh.magic[0] = PIF_MAGIC_NUM;
00306                 pfh.magic[1] = PIF_MAGIC_NUM;
00307                 sprintf(pfh.scalefactor, "1.0");
00308 
00309                 pfh.mode = PIF_FLOAT_INT;
00310                 sprintf(pfh.program, "EMAN %d/%02d/%02d", t->tm_mon, t->tm_mday, t->tm_year);
00311 
00312                 pfh.htype = 1;
00313                 pfh.nx = dict["nx"];
00314                 pfh.ny = dict["ny"];
00315                 pfh.nz = dict["nz"];
00316                 pfh.nimg += 1;
00317                 pfh.endian = (int) ByteOrder::is_host_big_endian();
00318                 fwrite(&pfh, sizeof(PifFileHeader), 1, pif_file);
00319         }
00320 
00321         PifImageHeader pih;
00322         memset(&pih, 0, sizeof(PifImageHeader));
00323         pih.nx = dict["nx"];
00324         pih.ny = dict["ny"];
00325         pih.nz = dict["nz"];
00326 
00327         pih.mode = PIF_FLOAT;
00328         pih.xlen = dict["apix_x"];
00329         pih.ylen = dict["apix_y"];
00330         pih.zlen = dict["apix_z"];
00331         pih.alpha = 90;
00332         pih.beta = 90;
00333         pih.gamma = 90;
00334         pih.mapc = 1;
00335         pih.mapr = 2;
00336         pih.maps = 3;
00337         pih.min = dict["minimum"];
00338         pih.max = dict["maximum"];
00339         pih.mean = dict["mean"];
00340         pih.sigma = dict["sigma"];
00341 
00342         pih.xorigin = dict["origin_x"];
00343         pih.yorigin = dict["origin_y"];
00344 
00345         sprintf(pih.time, "%d/%02d/%02d %d:%02d",
00346                         t->tm_mon, t->tm_mday, t->tm_year, t->tm_hour, t->tm_min);
00347         fwrite(&pih, sizeof(PifImageHeader), 1, pif_file);
00348 
00349         EXITFUNC;
00350         return 0;
00351 }
00352 
00353 int PifIO::read_data(float *data, int image_index, const Region *area, bool)
00354 {
00355         ENTERFUNC;
00356 
00357         check_read_access(image_index, data);
00358         fseek_to(image_index);
00359 
00360         int pih_sz = sizeof(PifImageHeader);
00361         PifImageHeader pih;
00362 
00363         if (fread(&pih, pih_sz, 1, pif_file) != 1) {
00364                 throw ImageReadException(filename, "PIF Image header");
00365         }
00366 
00367         if (area) {
00368                 check_region(area, FloatSize(pih.nx, pih.ny, pih.nz), is_new_file);
00369         }
00370 
00371         PifDataMode data_mode = static_cast < PifDataMode > (pih.mode);
00372         int num_layers = pih.nz;
00373 #if 0
00374         if (pfh.nz == pfh.nimg) {
00375                 num_layers = 1;
00376         }
00377 #endif
00378         // new way to read PIF data. The new way includes region reading.
00379         // If it is tested to be OK, remove the code in #if 0 ... #endif
00380         unsigned char * cdata = (unsigned char*)data;
00381         short *sdata = (short*) data;
00382 
00383         EMUtil::process_region_io(cdata, pif_file, READ_ONLY,
00384                                                           0, mode_size, pih.nx, pih.ny,
00385                                                           num_layers, area);
00386 
00387         int xlen = 0, ylen = 0, zlen = 0;
00388         EMUtil::get_region_dims(area, pih.nx, &xlen, pih.ny, &ylen, pih.nz, &zlen);
00389         size_t size = xlen * ylen * zlen;
00390 
00391         if(data_mode == PIF_FLOAT || data_mode == PIF_FLOAT_COMPLEX)
00392         {
00393                 become_host_endian< float >(data, size);
00394         }
00395         else {
00396                 if (mode_size == sizeof(short)) {
00397                         become_host_endian((short *) sdata, size);
00398                 }
00399                 else if (mode_size == sizeof(int)) {
00400                         become_host_endian((int *) data, size);
00401                 }
00402 
00403                 if (mode_size == sizeof(char)) {
00404                         for (size_t i = 0; i < size; i++) {
00405                                 size_t j = size - 1 - i;
00406                                 data[j] = (float)(cdata[j]) * real_scale_factor;
00407                         }
00408                 }
00409                 else if (mode_size == sizeof(short)) {
00410                         for (size_t i = 0; i < size; i++) {
00411                                 size_t j = size - 1 - i;
00412                                 data[j] = (float)(sdata[j]) * real_scale_factor;
00413                         }
00414                 }
00415                 else if (mode_size == sizeof(int)) {
00416                         for (size_t i = 0; i < size; i++) {
00417                                 size_t j = size - 1 - i;
00418                                 data[j] = (float) ((int *)data)[j] * real_scale_factor;
00419                         }
00420                 }
00421         }
00422 
00423         // end of new way for region reading
00424 
00425 #if 0
00426         int buf_size = pih.nx * mode_size;
00427         unsigned char *buf = new unsigned char[buf_size];
00428 
00429         for (int l = 0; l < num_layers; l++) {
00430                 int offset1 = l * pfh.nx * pfh.ny;
00431                 for (int j = 0; j < pfh.ny; j++) {
00432                         if (fread(buf, mode_size, pfh.nx, pif_file) != (unsigned int) pfh.nx) {
00433                                 if( buf )
00434                                 {
00435                                         delete[]buf;
00436                                         buf = 0;
00437                                 }
00438                                 throw ImageReadException(filename, "incomplete PIF read");
00439                         }
00440 
00441                         if (mode_size == sizeof(short)) {
00442                                 become_host_endian((short *) buf, pfh.nx);
00443                         }
00444                         else if (mode_size == sizeof(int)) {
00445                                 become_host_endian((int *) buf, pfh.nx);
00446                         }
00447 
00448                         int offset2 = offset1 + j * pfh.nx;
00449 
00450                         for (int k = 0; k < pfh.nx; k++) {
00451                                 float curr_data = 0;
00452 
00453                                 if (mode_size == sizeof(char)) {
00454                                         curr_data = (float) buf[k];
00455                                 }
00456                                 else if (mode_size == sizeof(short)) {
00457                                         curr_data = (float) ((short *) buf)[k];
00458                                 }
00459                                 else if (mode_size == sizeof(int)) {
00460                                         curr_data = (float) ((int *) buf)[k];
00461                                 }
00462                                 data[offset2 + k] = curr_data * real_scale_factor;
00463                         }
00464                 }
00465         }
00466         if( buf )
00467         {
00468                 delete[]buf;
00469                 buf = 0;
00470         }
00471 #endif
00472         EXITFUNC;
00473         return 0;
00474 }
00475 
00476 
00477 int PifIO::write_data(float *data, int image_index, const Region* area,
00478                                           EMUtil::EMDataType, bool)
00479 {
00480         ENTERFUNC;
00481 
00482         check_write_access(rw_mode, image_index, 0, data);
00483         fseek_to(image_index);
00484 
00485         int nx = pfh.nx;
00486         int ny = pfh.ny;
00487         int nz = pfh.nz;
00488 
00489         check_region(area, FloatSize(nx, ny, nz), is_new_file);
00490 
00491         // to write a region; if it works, please remove the code
00492         // in #if 0 ... #endif
00493         EMUtil::process_region_io(data, pif_file, WRITE_ONLY, 0,
00494                                                           mode_size, nx, ny, nz, area);
00495 
00496 #if 0
00497         size_t idx;
00498         int *buf = new int[nx];
00499         for (int i = 0; i < nz; i++) {
00500                 for (int j = 0; j < ny; j++) {
00501                         for (int k = 0; k < pfh.nx; k++) {
00502                                 idx = i * nx * ny + j * nx + k;
00503                                 buf[k] = (int) data[idx];
00504                         }
00505                         fwrite(buf, sizeof(int) * nx, 1, pif_file);
00506                 }
00507         }
00508         if( buf )
00509         {
00510                 delete[]buf;
00511                 buf = 0;
00512         }
00513 #endif
00514 
00515         EXITFUNC;
00516         return 0;
00517 }
00518 
00519 void PifIO::flush()
00520 {
00521         fflush(pif_file);
00522 }
00523 
00524 bool PifIO::is_complex_mode()
00525 {
00526         init();
00527         if (pfh.mode == PIF_SHORT_COMPLEX ||
00528                 pfh.mode == PIF_FLOAT_INT_COMPLEX ||
00529                 pfh.mode == PIF_FLOAT_COMPLEX || pfh.mode == PIF_SHORT_FLOAT_COMPLEX) {
00530                 return true;
00531         }
00532         return false;
00533 }
00534 
00535 bool PifIO::is_image_big_endian()
00536 {
00537         init();
00538         return is_big_endian;
00539 }
00540 
00541 int PifIO::get_nimg()
00542 {
00543         init();
00544         return pfh.nimg;
00545 }
00546 
00547 
00548 int PifIO::to_em_datatype(int p)
00549 {
00550         PifDataMode mode = static_cast < PifDataMode > (p);
00551         EMUtil::EMDataType e = EMUtil::EM_UNKNOWN;
00552 
00553         switch (mode) {
00554         case PIF_CHAR:
00555         case PIF_BOXED_DATA:
00556                 e = EMUtil::EM_CHAR;
00557                 break;
00558 
00559         case PIF_SHORT:
00560         case PIF_SHORT_FLOAT:
00561         case PIF_MAP_FLOAT_SHORT:
00562                 e = EMUtil::EM_SHORT;
00563                 break;
00564 
00565         case PIF_SHORT_COMPLEX:
00566         case PIF_SHORT_FLOAT_COMPLEX:
00567                 e = EMUtil::EM_SHORT_COMPLEX;
00568                 break;
00569 
00570         case PIF_FLOAT:
00571         case PIF_FLOAT_INT:
00572         case PIF_MAP_FLOAT_INT:
00573                 e = EMUtil::EM_FLOAT;
00574                 break;
00575         case PIF_FLOAT_COMPLEX:
00576         case PIF_FLOAT_INT_COMPLEX:
00577                 e = EMUtil::EM_FLOAT_COMPLEX;
00578                 break;
00579         case PIF_INVALID:
00580                 e = EMUtil::EM_UNKNOWN;
00581                 break;
00582         }
00583         return e;
00584 }
00585 
00586 int PifIO::to_pif_datatype(int e)
00587 {
00588         PifDataMode m = PIF_INVALID;
00589 
00590         switch (e) {
00591         case EMUtil::EM_CHAR:
00592                 m = PIF_BOXED_DATA;
00593                 break;
00594         case EMUtil::EM_SHORT:
00595                 m = PIF_SHORT;
00596                 break;
00597         case EMUtil::EM_SHORT_COMPLEX:
00598                 m = PIF_SHORT_COMPLEX;
00599                 break;
00600         case EMUtil::EM_FLOAT:
00601                 m = PIF_FLOAT_INT;
00602                 break;
00603         case EMUtil::EM_FLOAT_COMPLEX:
00604                 m = PIF_FLOAT_COMPLEX;
00605                 break;
00606         default:
00607                 LOGERR("unknown PIF mode: %d", e);
00608         }
00609 
00610         return m;
00611 }

Generated on Tue Jul 12 13:48:59 2011 for EMAN2 by  doxygen 1.3.9.1