00001
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036 #include <cstring>
00037 #include "xplorio.h"
00038 #include "util.h"
00039 #include "emassert.h"
00040 #include "portable_fileio.h"
00041 #include "geometry.h"
00042
00043 #ifdef WIN32
00044 #include <time.h>
00045 #include <windows.h>
00046 #define MAXPATHLEN (MAX_PATH*4)
00047 #else
00048 #include <sys/param.h>
00049 #endif
00050
00051
00052 using namespace EMAN;
00053
00054 const string XplorIO::SECTION_MODE = "ZYX";
00055 const int XplorIO::NFLOAT_PER_LINE = 6;
00056 const int XplorIO::INTEGER_SIZE = 8;
00057 const int XplorIO::FLOAT_SIZE = 12;
00058 const char * XplorIO::OUTFORMAT = "%12.5E";
00059
00060 XplorIO::XplorIO(const string & file, IOMode rw)
00061 : filename(file), rw_mode(rw), xplor_file(0), initialized(false)
00062 {
00063 is_big_endian = ByteOrder::is_host_big_endian();
00064 is_new_file = false;
00065 nlines_in_header = 0;
00066
00067 nx = 0;
00068 ny = 0;
00069 nz = 0;
00070
00071 apix_x = 0;
00072 apix_y = 0;
00073 apix_z = 0;
00074
00075 cell_alpha = 0;
00076 cell_beta = 0;
00077 cell_gama = 0;
00078 }
00079
00080 XplorIO::~XplorIO()
00081 {
00082 if (xplor_file) {
00083 fclose(xplor_file);
00084 xplor_file = 0;
00085 }
00086 }
00087
00088 void XplorIO::init()
00089 {
00090 if (initialized) {
00091 return;
00092 }
00093
00094 ENTERFUNC;
00095 initialized = true;
00096 xplor_file = sfopen(filename, rw_mode, &is_new_file);
00097
00098 if (!is_new_file) {
00099 char first_block[1024];
00100 fread(&first_block, sizeof(char), sizeof(first_block), xplor_file);
00101 if (!is_valid(&first_block)) {
00102 throw ImageReadException(filename, "invalid XPLOR");
00103 }
00104 portable_fseek(xplor_file, 0, SEEK_SET);
00105 char line[1024];
00106 int i = 1;
00107 int ntitle = 0;
00108
00109 int xmin = 0;
00110 int xmax = 0;
00111 int ymin = 0;
00112 int ymax = 0;
00113 int zmin = 0;
00114 int zmax = 0;
00115
00116 float cellx = 0;
00117 float celly = 0;
00118 float cellz = 0;
00119
00120 while(fgets(line, sizeof(line), xplor_file)) {
00121 line[strlen(line)-1] = '\0';
00122 if (i == 2) {
00123 ntitle = atoi(line);
00124 }
00125 else if (i == (ntitle+3)) {
00126 if (sscanf(line, "%8d%8d%8d%8d%8d%8d%8d%8d%8d", &nx, &xmin, &xmax,
00127 &ny, &ymin, &ymax, &nz, &zmin, &zmax) != 9) {
00128 throw ImageReadException(filename, "invalid XPLOR");
00129 }
00130 }
00131 else if (i == (ntitle+4)) {
00132 if(sscanf(line, "%f %f %f %f %f %f",
00133 &cellx, &celly, &cellz, &cell_alpha, &cell_beta, &cell_gama) != 6) {
00134 throw ImageReadException(filename, "invalid XPLOR");
00135 }
00136 }
00137 else if (i == (ntitle+5)) {
00138 break;
00139 }
00140
00141 i++;
00142 }
00143 nlines_in_header = i;
00144 apix_x = cellx / nx;
00145 apix_y = celly / ny;
00146 apix_z = cellz / nz;
00147 }
00148
00149 EXITFUNC;
00150 }
00151
00152 bool XplorIO::is_valid(const void *first_block)
00153 {
00154 ENTERFUNC;
00155 if (!first_block) {
00156 return false;
00157 }
00158 char *buf = (char *)(first_block);
00159 string line1 = Util::get_line_from_string(&buf);
00160 bool result = true;
00161
00162 if (line1.size() != 0) {
00163 result = false;
00164 }
00165 else {
00166 string line2 = Util::get_line_from_string(&buf);
00167 int ntitle = 0;
00168
00169 if ((int)line2.size() != INTEGER_SIZE) {
00170 result = false;
00171 }
00172 else {
00173 ntitle = atoi(line2.c_str());
00174 if (ntitle < 0 || ntitle > 50) {
00175 result = false;
00176 }
00177
00178 else {
00179 for (int i = 0; i < ntitle+2; i++) {
00180 Util::get_line_from_string(&buf);
00181 }
00182
00183 string modeline = Util::get_line_from_string(&buf);
00184 if (modeline != SECTION_MODE) {
00185 result = false;
00186 }
00187 }
00188 }
00189 }
00190
00191 EXITFUNC;
00192 return result;
00193 }
00194
00195 int XplorIO::read_header(Dict &dict, int image_index, const Region *area, bool)
00196 {
00197 ENTERFUNC;
00198
00199
00200 if(image_index == -1) {
00201 image_index = 0;
00202 }
00203 if(image_index != 0) {
00204 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().");
00205 }
00206
00207 init();
00208 check_region(area, FloatSize(nx, ny, nz), is_new_file);
00209
00210 int xlen = 0, ylen = 0, zlen = 0;
00211 EMUtil::get_region_dims(area, nx, &xlen, ny, &ylen, nz, &zlen);
00212
00213 dict["nx"] = xlen;
00214 dict["ny"] = ylen;
00215 dict["nz"] = zlen;
00216
00217 dict["apix_x"] = apix_x;
00218 dict["apix_y"] = apix_y;
00219 dict["apix_z"] = apix_z;
00220
00221 dict["XPLOR.alpha"] = cell_alpha;
00222 dict["XPLOR.beta"] = cell_beta;
00223 dict["XPLOR.gama"] = cell_gama;
00224
00225 EXITFUNC;
00226 return 0;
00227 }
00228
00229 int XplorIO::write_header(const Dict & dict, int image_index, const Region* area,
00230 EMUtil::EMDataType, bool)
00231 {
00232 ENTERFUNC;
00233
00234 if(image_index == -1) {
00235 image_index = 0;
00236 }
00237 if(image_index != 0) {
00238 throw ImageWriteException(filename, "XPLOR file does not support stack.");
00239 }
00240 check_write_access(rw_mode, image_index);
00241 if (area) {
00242 check_region(area, FloatSize(nx, ny, nz), is_new_file);
00243 EXITFUNC;
00244 return 0;
00245 }
00246
00247 nx = dict["nx"];
00248 ny = dict["ny"];
00249 nz = dict["nz"];
00250 apix_x = dict["apix_x"];
00251 apix_y = dict["apix_y"];
00252 apix_z = dict["apix_z"];
00253
00254 nlines_in_header = 0;
00255 time_t t0 = time(0);
00256 struct tm *t = localtime(&t0);
00257 rewind(xplor_file);
00258
00259 fprintf(xplor_file, "\n");
00260 fprintf(xplor_file, "%8d\n", 1);
00261 fprintf(xplor_file, "\"%s\" written by EMAN at %s", filename.c_str(), asctime(t));
00262
00263 int z0 = -nz / 2;
00264 int z1 = (nz - 1) / 2;
00265
00266 if (2 * nz - 1 == nx && 2 * nz - 1 == ny) {
00267 z0 = 0;
00268 z1 = nz - 1;
00269 }
00270
00271 fprintf(xplor_file, "%8d%8d%8d%8d%8d%8d%8d%8d%8d\n",
00272 nx, -nx / 2, nx % 2 ? nx / 2 : nx / 2 - 1, ny, -ny / 2,
00273 ny % 2 ? ny / 2 : ny / 2 - 1, nz, z0, z1);
00274
00275 char fformat[256];
00276 sprintf(fformat, "%s%s%s%s%s%s\n",
00277 OUTFORMAT, OUTFORMAT,OUTFORMAT, OUTFORMAT, OUTFORMAT,OUTFORMAT);
00278
00279 fprintf(xplor_file, fformat,
00280 nx * apix_x, ny * apix_y, nz * apix_z, 90.0, 90.0, 90.0);
00281 fprintf(xplor_file, "ZYX\n");
00282 nlines_in_header = 5;
00283 flush();
00284
00285 EXITFUNC;
00286 return 0;
00287 }
00288
00289 int XplorIO::read_data(float *data, int image_index, const Region *area, bool)
00290 {
00291 ENTERFUNC;
00292
00293
00294 image_index = 0;
00295 check_read_access(image_index);
00296 FloatSize max_size = FloatSize(nx, ny, nz);
00297 check_region(area, max_size, is_new_file);
00298
00299
00300
00301
00302
00303
00304 if (area != 0 && !area->is_region_in_box(max_size)) {
00305 char desc[1024];
00306 sprintf(desc, "Region box %s is outside image area (%d,%d,%d)",
00307 area->get_string().c_str(), (int)max_size[0],
00308 (int)max_size[1], (int)max_size[2]);
00309 throw ImageReadException("", desc);
00310 }
00311
00312 Assert(nlines_in_header > 0);
00313 rewind(xplor_file);
00314 EMUtil::jump_lines(xplor_file, nlines_in_header);
00315
00316 EMUtil::process_ascii_region_io(data, xplor_file, ImageIO::READ_ONLY, image_index,
00317 FLOAT_SIZE, nx, ny, nz, area, true,
00318 NFLOAT_PER_LINE, OUTFORMAT);
00319
00320 EXITFUNC;
00321 return 0;
00322 }
00323
00324
00325 #if 0
00326 int XplorIO::read_data(float *data, int, const Region *, bool)
00327 {
00328 ENTERFUNC;
00329 int step = NFLOAT_PER_LINE;
00330 char line[1024];
00331 int nxy = nx * ny;
00332 int nlines = nxy / step;
00333 int nleft = nxy - nlines * step;
00334
00335 for (int k = 0; k < nz; k++) {
00336 fgets(line, sizeof(line), xplor_file);
00337 int kk = 0;
00338 sscanf(line, "%d", &kk);
00339 if (kk != (k+1)) {
00340 LOGERR("section index = %d. It should be %d\n", kk, (k+1));
00341 }
00342
00343 int k2 = k * nxy;
00344
00345 for (int i = 0; i < nlines; i++) {
00346 fgets(line, sizeof(line), xplor_file);
00347 int i2 = k2 + i * step;
00348 sscanf(line, "%f %f %f %f %f %f",
00349 &data[i2], &data[i2+1], &data[i2+2],
00350 &data[i2+3], &data[i2+4], &data[i2+5]);
00351 }
00352
00353 if (nleft > 0) {
00354 int i2 = k2 + nlines * step;
00355 fgets(line, sizeof(line), xplor_file);
00356 char *pline = line;
00357 for (int j = 0; j < nleft; j++) {
00358 sscanf(pline, "%f", &data[i2+j]);
00359 pline += FLOAT_SIZE;
00360 }
00361 }
00362 }
00363
00364
00365 EXITFUNC;
00366 return 0;
00367 }
00368 #endif
00369
00370
00371 int XplorIO::write_data(float *data, int image_index, const Region* area,
00372 EMUtil::EMDataType, bool)
00373 {
00374
00375 ENTERFUNC;
00376
00377 image_index = 0;
00378 check_write_access(rw_mode, image_index, 1, data);
00379 check_region(area, FloatSize(nx,ny,nz), is_new_file);
00380
00381 if (!is_new_file) {
00382 rewind(xplor_file);
00383 EMUtil::jump_lines(xplor_file, nlines_in_header);
00384 }
00385
00386 int nsecs = nx * ny;
00387 int step = NFLOAT_PER_LINE;
00388
00389 if (!area) {
00390 for (int k = 0; k < nz; ++k) {
00391 fprintf(xplor_file, "%8d\n", (k+1));
00392
00393 for (int i = 0; i < nsecs - step; i += step) {
00394 for (int j = 0; j < step; j++) {
00395 fprintf(xplor_file, OUTFORMAT, data[k * nsecs + i + j]);
00396 }
00397 fprintf(xplor_file, "\n");
00398 }
00399
00400 for (int l = (nsecs - 1) / step * step; l < nsecs; l++) {
00401 fprintf(xplor_file, OUTFORMAT, data[k * nsecs + l]);
00402 }
00403
00404 fprintf(xplor_file, "\n");
00405 }
00406
00407
00408
00409 }
00410 else {
00411 EMUtil::process_ascii_region_io(data, xplor_file, ImageIO::WRITE_ONLY,
00412 image_index, FLOAT_SIZE, nx, ny, nz,
00413 area, true, NFLOAT_PER_LINE, OUTFORMAT);
00414
00415 }
00416
00417 EXITFUNC;
00418
00419 return 0;
00420 }
00421
00422 void XplorIO::flush()
00423 {
00424 fflush(xplor_file);
00425 }
00426
00427 bool XplorIO::is_complex_mode()
00428 {
00429 return false;
00430 }
00431
00432 bool XplorIO::is_image_big_endian()
00433 {
00434 init();
00435 return is_big_endian;
00436 }