ecat63r.c
Go to the documentation of this file.
00001 /******************************************************************************
00002 
00003   ecat63r.c  (c) 2003,2004 by Turku PET Centre
00004 
00005   Procedures for reading ECAT 6.3 format.
00006 
00007   Assumptions:
00008   1. Data_type in headers specifies, whether ints, long ints and floats in
00009      header and matrix data are in VAX format (1, 2, 3 and 4) or in IEEE
00010      (i386, SUN) format.
00011   2. Data is automatically converted to little or big endian when read,
00012      according to the current platform.
00013   3. Data is automatically converted out from the VAX format when read,
00014      but header data_type remains as original.
00015 
00016   Version: Contents were previously in ecat63.c
00017   2003-07-21 Vesa Oikonen
00018   2003-09-05 VO
00019     Added function ecat63ReadImageMatrix() and ecat63ReadScanMatrix().
00020   2003-09-08 VO
00021     Added function ecat63pxlbytes().
00022   2004-06-21 VO
00023     Larger range allowed for calibration factor in ecat63ReadMainheader().
00024   2004-09-20 VO
00025     Doxygen style comments.
00026 
00027 ******************************************************************************/
00028 #include <stdio.h>
00029 #include <stdlib.h>
00030 #include <math.h>
00031 #include <ctype.h>
00032 #include <string.h>
00033 #include <unistd.h>
00034 #include <time.h>
00035 /*****************************************************************************/
00036 #include "swap.h"
00037 #include "include/ecat63.h"
00038 /*****************************************************************************/
00039 
00040 /*****************************************************************************/
00050 int ecat63ReadMainheader(FILE *fp, ECAT63_mainheader *h) {
00051   unsigned char buf[MatBLKSIZE];
00052   int i;
00053   int little; /* 1 if current platform is little endian (i386), else 0 */
00054   int vaxdata=1; /* 1 if data is in VAX format, else 0 */
00055 
00056   if(ECAT63_TEST) printf("ecat63ReadMainheader()\n");
00057   if(fp==NULL || h==NULL) return(1);
00058   little=little_endian();
00059   /* Seek the first block */
00060   fseek(fp, 0, SEEK_SET); if(ftell(fp)!=0) return(1);
00061   /* Read the subheader block */
00062   if(fread(buf, MatBLKSIZE, 1, fp)<1) return(2);
00063 
00064   /* Copy char data to header structure */
00065   memcpy(h->ecat_format, buf+0, 14); memcpy(h->fill1, buf+14, 14);
00066   memcpy(h->original_file_name, buf+28, 20); memcpy(h->node_id, buf+56, 10);
00067   memcpy(h->isotope_code, buf+78, 8); memcpy(h->radiopharmaceutical, buf+90, 32);
00068   memcpy(h->study_name, buf+162, 12); memcpy(h->patient_id, buf+174, 16);
00069   memcpy(h->patient_name, buf+190, 32); h->patient_sex=buf[222];
00070   memcpy(h->patient_age, buf+223, 10); memcpy(h->patient_height, buf+233, 10);
00071   memcpy(h->patient_weight, buf+243, 10); h->patient_dexterity=buf[253];
00072   memcpy(h->physician_name, buf+254, 32); memcpy(h->operator_name, buf+286, 32);
00073   memcpy(h->study_description, buf+318, 32); memcpy(h->facility_name, buf+356, 20);
00074   memcpy(h->user_process_code, buf+462, 10);
00075 
00076   /* Copy short ints */
00077   /* in big endian platform, change byte order temporarily */
00078   if(!little) swabip(buf, MatBLKSIZE);
00079   memcpy(&h->data_type, buf+50, 2); if(h->data_type>4) vaxdata=0;
00080   /*printf("main_header.data_type=%d\n", h->data_type);*/
00081   memcpy(&h->sw_version, buf+48, 2);
00082   memcpy(&h->system_type, buf+52, 2); memcpy(&h->file_type, buf+54, 2);
00083   memcpy(&h->scan_start_day, buf+66, 2); memcpy(&h->scan_start_month, buf+68, 2);
00084   memcpy(&h->scan_start_year, buf+70, 2); memcpy(&h->scan_start_hour, buf+72, 2);
00085   memcpy(&h->scan_start_minute, buf+74, 2); memcpy(&h->scan_start_second, buf+76, 2);
00086   memcpy(&h->rot_source_speed, buf+134, 2); memcpy(&h->wobble_speed, buf+136, 2);
00087   memcpy(&h->transm_source_type, buf+138, 2); memcpy(&h->transaxial_samp_mode, buf+148, 2);
00088   memcpy(&h->coin_samp_mode, buf+150, 2); memcpy(&h->axial_samp_mode, buf+152, 2);
00089   memcpy(&h->calibration_units, buf+158, 2); memcpy(&h->compression_code, buf+160, 2);
00090   memcpy(&h->acquisition_type, buf+350, 2); memcpy(&h->bed_type, buf+352, 2);
00091   memcpy(&h->septa_type, buf+354, 2); memcpy(&h->num_planes, buf+376, 2);
00092   memcpy(&h->num_frames, buf+378, 2); memcpy(&h->num_gates, buf+380, 2);
00093   memcpy(&h->num_bed_pos, buf+382, 2); memcpy(&h->lwr_sctr_thres, buf+452, 2);
00094   memcpy(&h->lwr_true_thres, buf+454, 2); memcpy(&h->upr_true_thres, buf+456, 2);
00095   memcpy(h->fill2, buf+472, 40);
00096   /* Change back the byte order */
00097   if(!little) swabip(buf, MatBLKSIZE);
00098 
00099   /* Copy floats */
00100   h->isotope_halflife=ecat63rFloat(buf+86, vaxdata, little);
00101   h->gantry_tilt=ecat63rFloat(buf+122, vaxdata, little);
00102   h->gantry_rotation=ecat63rFloat(buf+126, vaxdata, little);
00103   h->bed_elevation=ecat63rFloat(buf+130, vaxdata, little);
00104   h->axial_fov=ecat63rFloat(buf+140, vaxdata, little);
00105   h->transaxial_fov=ecat63rFloat(buf+144, vaxdata, little);
00106   h->calibration_factor=ecat63rFloat(buf+154, vaxdata, little);
00107   h->init_bed_position=ecat63rFloat(buf+384, vaxdata, little);
00108   for(i=0; i<15; i++) h->bed_offset[i]=ecat63rFloat(buf+388+i*4, vaxdata, little);
00109   h->plane_separation=ecat63rFloat(buf+448, vaxdata, little);
00110   h->collimator=ecat63rFloat(buf+458, vaxdata, little);
00111 
00112   /* Check file format and platform */
00113   if(ECAT63_TEST) printf("ecat_format='%.14s'\n", h->ecat_format);
00114   /* if format is not specified, ECAT63 is assumed */
00115   if(h->ecat_format[0]==(char)0) {
00116     strcpy(h->ecat_format, "ECAT63");
00117   }
00118   /* only ECAT63 format is approved here */
00119   if(ECAT63_TEST) printf("ecat_format='%.14s'\n", h->ecat_format);
00120   if(strncmp(h->ecat_format, "ECAT63", 6)!=0) return(3);
00121 
00122   /* Check that most important contents are ok */
00123   if(h->data_type<BYTE_TYPE || h->data_type>SUN_I4) {
00124     if(ECAT63_TEST) printf("Invalid data types; probable conversion error.\n");
00125     return(5);
00126   }
00127   if(h->calibration_factor<0.0 || h->calibration_factor>1.0e12) {
00128     if(ECAT63_TEST) printf("Invalid calibration factor; possible conversion error.\n");
00129     return(6);
00130   }
00131   if(h->file_type!=RAW_DATA && h->file_type!=IMAGE_DATA &&
00132      h->file_type!=ATTN_DATA && h->file_type!=NORM_DATA) {
00133     if(ECAT63_TEST) printf("Invalid file types; probable conversion error.\n");
00134     return(7);
00135   }
00136 
00137   return(0);
00138 }
00139 /*****************************************************************************/
00140 
00141 /*****************************************************************************/
00152 int ecat63ReadImageheader(FILE *fp, int blk, ECAT63_imageheader *h) {
00153   unsigned char buf[MatBLKSIZE];
00154   int i;
00155   int little; /* 1 if current platform is little endian (i386), else 0 */
00156   int vaxdata=1; /* 1 if data is in VAX format, else 0 */
00157 
00158   if(ECAT63_TEST) printf("ecat63ReadImageheader(fp, %d ih)\n", blk);
00159   if(fp==NULL || blk<2 || h==NULL) return(1);
00160   little=little_endian();
00161   /* Seek the subheader block */
00162   fseek(fp, (blk-1)*MatBLKSIZE, SEEK_SET); if(ftell(fp)!=(blk-1)*MatBLKSIZE) return(2);
00163   /* Read the subheader block */
00164   if(fread(buf, MatBLKSIZE, 1, fp)<1) return(3);
00165 
00166   /* Copy char data to header structure */
00167   memcpy(h->fill1, buf+0, 126); memcpy(h->annotation, buf+420, 40);
00168 
00169   /* Copy short ints */
00170   /* in big endian platform, change byte order temporarily */
00171   if(!little) swabip(buf, MatBLKSIZE);
00172   memcpy(&h->data_type, buf+126, 2); if(h->data_type>4) vaxdata=0;
00173   /*printf("data_type=%d\n", h->data_type);*/
00174   memcpy(&h->num_dimensions, buf+128, 2);
00175   memcpy(&h->dimension_1, buf+132, 2); memcpy(&h->dimension_2, buf+134, 2);
00176   memcpy(&h->image_min, buf+176, 2); memcpy(&h->image_max, buf+178, 2);
00177   memcpy(&h->slice_location, buf+200, 2); memcpy(&h->recon_start_hour, buf+202, 2);
00178   memcpy(&h->recon_start_min, buf+204, 2); memcpy(&h->recon_start_sec, buf+206, 2);
00179   memcpy(&h->filter_code, buf+236, 2); memcpy(&h->processing_code, buf+376, 2);
00180   memcpy(&h->quant_units, buf+380, 2); memcpy(&h->recon_start_day, buf+382, 2);
00181   memcpy(&h->recon_start_month, buf+384, 2); memcpy(&h->recon_start_year, buf+386, 2);
00182   memcpy(h->fill2, buf+460, 52);
00183   /* Change back the byte order */
00184   if(!little) swabip(buf, MatBLKSIZE);
00185 
00186   /* Copy ints */
00187   h->frame_duration=ecat63rInt(buf+192, vaxdata, little);
00188   h->frame_start_time=ecat63rInt(buf+196, vaxdata, little);
00189   h->recon_duration=ecat63rInt(buf+208, vaxdata, little);
00190   h->scan_matrix_num=ecat63rInt(buf+238, vaxdata, little);
00191   h->norm_matrix_num=ecat63rInt(buf+242, vaxdata, little);
00192   h->atten_cor_mat_num=ecat63rInt(buf+246, vaxdata, little);
00193 
00194   /* Copy floats */
00195   h->x_origin=ecat63rFloat(buf+160, vaxdata, little);
00196   h->y_origin=ecat63rFloat(buf+164, vaxdata, little);
00197   h->recon_scale=ecat63rFloat(buf+168, vaxdata, little);
00198   h->quant_scale=ecat63rFloat(buf+172, vaxdata, little);
00199   h->pixel_size=ecat63rFloat(buf+184, vaxdata, little);
00200   h->slice_width=ecat63rFloat(buf+188, vaxdata, little);
00201   h->image_rotation=ecat63rFloat(buf+296, vaxdata, little);
00202   h->plane_eff_corr_fctr=ecat63rFloat(buf+300, vaxdata, little);
00203   h->decay_corr_fctr=ecat63rFloat(buf+304, vaxdata, little);
00204   h->loss_corr_fctr=ecat63rFloat(buf+308, vaxdata, little);
00205   h->intrinsic_tilt=ecat63rFloat(buf+312, vaxdata, little);
00206   h->ecat_calibration_fctr=ecat63rFloat(buf+388, vaxdata, little);
00207   h->well_counter_cal_fctr=ecat63rFloat(buf+392, vaxdata, little);
00208   for(i=0; i<6; i++) h->filter_params[i]=ecat63rFloat(buf+396+i*4, vaxdata, little);
00209 
00210   /* Check that most important contents are ok */
00211   if(h->data_type<BYTE_TYPE || h->data_type>SUN_I4) {
00212     if(ECAT63_TEST) printf("Invalid data types; probable conversion error.\n");
00213     return(4);
00214   }
00215   if(h->ecat_calibration_fctr<0.0 || h->ecat_calibration_fctr>1.0e10) {
00216     if(ECAT63_TEST) printf("Invalid calibration factor; probable conversion error.\n");
00217     return(5);
00218   }
00219   if(h->frame_duration<0.0 || h->frame_duration>1.0e12) {
00220     if(ECAT63_TEST) printf("Invalid frame duration; probable conversion error.\n");
00221     return(6);
00222   }
00223 
00224   return(0);
00225 }
00226 /*****************************************************************************/
00227 
00228 /*****************************************************************************/
00238 int ecat63ReadAttnheader(FILE *fp, int blk, ECAT63_attnheader *h) {
00239   unsigned char buf[MatBLKSIZE];
00240   int little; /* 1 if current platform is little endian (i386), else 0 */
00241   int vaxdata=1; /* 1 if data is in VAX format, else 0 */
00242 
00243   if(ECAT63_TEST) printf("ecat63ReadAttnheader(fp, %d, ah)\n", blk);
00244   if(fp==NULL || blk<2 || h==NULL) return(1);
00245   little=little_endian();
00246   /* Seek the subheader block */
00247   fseek(fp, (blk-1)*MatBLKSIZE, SEEK_SET); if(ftell(fp)!=(blk-1)*MatBLKSIZE) return(2);
00248   /* Read the subheader block */
00249   if(fread(buf, MatBLKSIZE, 1, fp)<1) return(3);
00250 
00251   /* Copy short ints */
00252   /* in big endian platform, change byte order temporarily */
00253   if(!little) swabip(buf, MatBLKSIZE);
00254   memcpy(&h->data_type, buf+126, 2); if(h->data_type>4) vaxdata=0;
00255   /*printf("data_type=%d\n", h->data_type);*/
00256   memcpy(&h->attenuation_type, buf+128, 2);
00257   memcpy(&h->dimension_1, buf+132, 2); memcpy(&h->dimension_2, buf+134, 2);
00258   /* Change back the byte order */
00259   if(!little) swabip(buf, MatBLKSIZE);
00260 
00261   /* Copy floats */
00262   h->scale_factor=ecat63rFloat(buf+182, vaxdata, little);
00263   h->x_origin=ecat63rFloat(buf+186, vaxdata, little);
00264   h->y_origin=ecat63rFloat(buf+190, vaxdata, little);
00265   h->x_radius=ecat63rFloat(buf+194, vaxdata, little);
00266   h->y_radius=ecat63rFloat(buf+198, vaxdata, little);
00267   h->tilt_angle=ecat63rFloat(buf+202, vaxdata, little);
00268   h->attenuation_coeff=ecat63rFloat(buf+206, vaxdata, little);
00269   h->sample_distance=ecat63rFloat(buf+210, vaxdata, little);
00270 
00271   /* Check that most important contents are ok */
00272   if(h->data_type<BYTE_TYPE || h->data_type>SUN_I4) {
00273     if(ECAT63_TEST) printf("Invalid data types; probable conversion error.\n");
00274     return(4);
00275   }
00276   if(h->scale_factor<=0.0 || h->scale_factor>1.0e8) {
00277     if(ECAT63_TEST) printf("Invalid scale factor; probable conversion error.\n");
00278     return(5);
00279   }
00280 
00281   return(0);
00282 }
00283 /*****************************************************************************/
00284 
00285 /*****************************************************************************/
00296 int ecat63ReadScanheader(FILE *fp, int blk, ECAT63_scanheader *h) {
00297   unsigned char buf[MatBLKSIZE];
00298   int i;
00299   int little; /* 1 if current platform is little endian (i386), else 0 */
00300   int vaxdata=1; /* 1 if data is in VAX format, else 0 */
00301 
00302   if(ECAT63_TEST) printf("ecat63ReadScanheader(fp, %d, sh)\n", blk);
00303   if(fp==NULL || blk<2 || h==NULL) return(1);
00304   little=little_endian();
00305   /* Seek the subheader block */
00306   fseek(fp, (blk-1)*MatBLKSIZE, SEEK_SET); if(ftell(fp)!=(blk-1)*MatBLKSIZE) return(2);
00307   /* Read the subheader block */
00308   if(fread(buf, MatBLKSIZE, 1, fp)<1) return(3);
00309 
00310   /* Copy char data to header structure */
00311   memcpy(h->fill1, buf+0, 126);
00312 
00313   /* Copy short ints */
00314   /* in big endian platform, change byte order temporarily */
00315   if(!little) swabip(buf, MatBLKSIZE);
00316   memcpy(&h->data_type, buf+126, 2); if(h->data_type>4) vaxdata=0;
00317   /*printf("data_type=%d\n", h->data_type);*/
00318   memcpy(&h->dimension_1, buf+132, 2); memcpy(&h->dimension_2, buf+134, 2);
00319   memcpy(&h->smoothing, buf+136, 2); memcpy(&h->processing_code, buf+138, 2);
00320   memcpy(&h->frame_duration_sec, buf+170, 2);
00321   memcpy(&h->scan_min, buf+192, 2); memcpy(&h->scan_max, buf+194, 2);
00322   memcpy(h->fill2, buf+468, 44);
00323   /* Change back the byte order */
00324   if(!little) swabip(buf, MatBLKSIZE);
00325 
00326   /* Copy ints */
00327   h->gate_duration=ecat63rInt(buf+172, vaxdata, little);
00328   h->r_wave_offset=ecat63rInt(buf+176, vaxdata, little);
00329   h->prompts=ecat63rInt(buf+196, vaxdata, little);
00330   h->delayed=ecat63rInt(buf+200, vaxdata, little);
00331   h->multiples=ecat63rInt(buf+204, vaxdata, little);
00332   h->net_trues=ecat63rInt(buf+208, vaxdata, little);
00333   h->total_coin_rate=ecat63rInt(buf+452, vaxdata, little);
00334   h->frame_start_time=ecat63rInt(buf+456, vaxdata, little);
00335   h->frame_duration=ecat63rInt(buf+460, vaxdata, little);
00336 
00337   /* Copy floats */
00338   h->sample_distance=ecat63rFloat(buf+146, vaxdata, little);
00339   h->isotope_halflife=ecat63rFloat(buf+166, vaxdata, little);
00340   h->scale_factor=ecat63rFloat(buf+182, vaxdata, little);
00341   for(i=0; i<16; i++) h->cor_singles[i]=ecat63rFloat(buf+316+i*4, vaxdata, little);
00342   for(i=0; i<16; i++) h->uncor_singles[i]=ecat63rFloat(buf+380+i*4, vaxdata, little);
00343   h->tot_avg_cor=ecat63rFloat(buf+444, vaxdata, little);
00344   h->tot_avg_uncor=ecat63rFloat(buf+448, vaxdata, little);
00345   h->loss_correction_fctr=ecat63rFloat(buf+464, vaxdata, little);
00346 
00347   /* Check that most important contents are ok */
00348   if(h->data_type<BYTE_TYPE || h->data_type>SUN_I4) {
00349     if(ECAT63_TEST) printf("Invalid data types; probable conversion error.\n");
00350     return(4);
00351   }
00352   if(h->scale_factor<=0.0 || h->scale_factor>1.0e8) {
00353     if(ECAT63_TEST) printf("Invalid scale factor; probable conversion error.\n");
00354     return(5);
00355   }
00356   if(h->frame_duration<0.0 || h->frame_duration>1.0e12) {
00357     if(ECAT63_TEST) printf("Invalid frame duration; probable conversion error.\n");
00358     return(6);
00359   }
00360 
00361   return(0);
00362 }
00363 /*****************************************************************************/
00364 
00365 /*****************************************************************************/
00375 int ecat63ReadNormheader(FILE *fp, int blk, ECAT63_normheader *h) {
00376   unsigned char buf[MatBLKSIZE];
00377   int little; /* 1 if current platform is little endian (i386), else 0 */
00378   int vaxdata=1; /* 1 if data is in VAX format, else 0 */
00379 
00380   if(ECAT63_TEST) printf("ecat63ReadNormheader(fp, %d, nh)\n", blk);
00381   if(fp==NULL || blk<2 || h==NULL) return(1);
00382   little=little_endian();
00383   /* Seek the subheader block */
00384   fseek(fp, (blk-1)*MatBLKSIZE, SEEK_SET); if(ftell(fp)!=(blk-1)*MatBLKSIZE) return(2);
00385   /* Read the subheader block */
00386   if(fread(buf, MatBLKSIZE, 1, fp)<1) return(3);
00387 
00388   /* Copy short ints */
00389   /* in big endian platform, change byte order temporarily */
00390   if(!little) swabip(buf, MatBLKSIZE);
00391   memcpy(&h->data_type, buf+126, 2); if(h->data_type>4) vaxdata=0;
00392   /*printf("data_type=%d\n", h->data_type);*/
00393   memcpy(&h->dimension_1, buf+132, 2); memcpy(&h->dimension_2, buf+134, 2);
00394   memcpy(&h->norm_hour, buf+186, 2); memcpy(&h->norm_minute, buf+188, 2);
00395   memcpy(&h->norm_second, buf+190, 2); memcpy(&h->norm_day, buf+192, 2);
00396   memcpy(&h->norm_month, buf+194, 2); memcpy(&h->norm_year, buf+196, 2);
00397   /* Change back the byte order */
00398   if(!little) swabip(buf, MatBLKSIZE);
00399 
00400   /* Copy floats */
00401   h->scale_factor=ecat63rFloat(buf+182, vaxdata, little);
00402   h->fov_source_width=ecat63rFloat(buf+198, vaxdata, little);
00403 
00404   /* Check that most important contents are ok */
00405   if(h->data_type<BYTE_TYPE && h->data_type>SUN_I4) {
00406     if(ECAT63_TEST) printf("Invalid data types; probable conversion error.\n");
00407     return(4);
00408   }
00409   if(h->scale_factor<=0.0 || h->scale_factor>1.0e8) {
00410     if(ECAT63_TEST) printf("Invalid scale factor; probable conversion error.\n");
00411     return(5);
00412   }
00413 
00414   return(0);
00415 }
00416 /*****************************************************************************/
00417 
00418 /*****************************************************************************/
00432 int ecat63ReadMatdata(FILE *fp, int strtblk, int blkNr, char *data, int dtype) {
00433   int i, n, little, err=0;
00434   char *cptr;
00435   float f;
00436 
00437 
00438   if(ECAT63_TEST) printf("ecat63ReadMatdata(fp, %d, %d, data, %d)\n", strtblk, blkNr, dtype);
00439   /* Check the arguments */
00440   if(blkNr<=0 || strtblk<1 || data==NULL) return(1);
00441   /* Seek the first data block */
00442   fseek(fp, (strtblk-1)*MatBLKSIZE, SEEK_SET);
00443   if(ftell(fp)!=(strtblk-1)*MatBLKSIZE) return(9);
00444   /* Read the data blocks */
00445   if(fread(data, MatBLKSIZE, blkNr, fp) < blkNr) return(2);
00446   /* Translate data if necessary */
00447   little=little_endian();
00448   switch(dtype) {
00449     case BYTE_TYPE: /* byte format...no translation necessary */
00450       break;
00451     case VAX_I2:    /* byte conversion necessary on big endian platform */
00452       if(!little) {cptr=data; swabip(cptr, blkNr*MatBLKSIZE);}
00453       break;
00454     case VAX_I4:
00455       for(i=0, cptr=data; i<blkNr*MatBLKSIZE; i+=4, cptr+=4) {
00456         n=ecat63rInt(cptr, 1, little); memcpy(cptr, &n, 4);
00457       }
00458       break;
00459     case VAX_R4:
00460       for(i=0, cptr=data; i<blkNr*MatBLKSIZE; i+=4, cptr+=4) {
00461         f=ecat63rFloat(cptr, 1, little); memcpy(cptr, &f, 4);
00462       }
00463       break;
00464     case IEEE_R4:   /* IEEE float ; byte conversion necessary on big end platforms */
00465     case SUN_I4:    /* SUN int ; byte conversion necessary on big end platforms */
00466       if(!little) swawbip(data, blkNr*MatBLKSIZE);
00467       break;
00468     case SUN_I2:    /* SUN short ; byte conversion necessary on big end platforms */
00469       if(!little) swabip(data, blkNr*MatBLKSIZE);
00470       break;
00471     default:  /* if something else, for now think it as an error */
00472       err=2;
00473       break;
00474   }
00475   return(err);
00476 }
00477 /*****************************************************************************/
00478 
00479 /*****************************************************************************/
00494 int ecat63ReadImageMatrix(FILE *fp, int first_block, int last_block, ECAT63_imageheader *h, float **fdata) {
00495   int i, ret, blockNr, pxlNr;
00496   char *mdata, *mptr;
00497   float *_fdata, *fptr;
00498   short int *sptr;
00499   int *iptr;
00500   
00501   
00502   if(ECAT63_TEST) printf("ecat63ReadImageMatrix(fp, %d, %d, hdr, fdata)\n",
00503     first_block, last_block);
00504   if(fp==NULL || first_block<=MatFirstDirBlk || h==NULL) {
00505     sprintf(ecat63errmsg, "invalid function parameter.\n");
00506     return(1);
00507   }
00508   *fdata=(float*)NULL;
00509   
00510   /* Read subheader */
00511   ret=ecat63ReadImageheader(fp, first_block, h);
00512   if(ret) {
00513     sprintf(ecat63errmsg, "cannot read subheader (%d).\n", ret);
00514     return(5);
00515   }
00516   if(ECAT63_TEST>4) ecat63PrintImageheader(h, stdout);
00517   pxlNr=h->dimension_1*h->dimension_2;
00518   if(pxlNr<=0) {
00519     sprintf(ecat63errmsg, "invalid matrix dimension.\n");
00520     return(6);  
00521   }
00522   
00523   /* Read matrix data */
00524   blockNr=last_block-first_block; if(blockNr<1) return(0);
00525   mdata=(char*)malloc(blockNr*MatBLKSIZE);
00526   if(mdata==NULL) {
00527     sprintf(ecat63errmsg, "cannot allocate memory.\n");
00528     return(8);  
00529   }
00530   mptr=mdata;
00531   ret=ecat63ReadMatdata(fp, first_block+1, blockNr, mptr, h->data_type);
00532   if(ret || mdata==NULL) {
00533     sprintf(ecat63errmsg, "cannot read matrix data (%d).\n", ret);
00534     free(mdata); return(9);
00535   }
00536   
00537   /* Allocate memory for float data */
00538   _fdata=(float*)malloc(pxlNr*sizeof(float));
00539   if(_fdata==NULL) {
00540     sprintf(ecat63errmsg, "cannot allocate memory.\n");
00541     free(mdata); return(11);  
00542   }
00543 
00544   /* Convert matrix data to floats */
00545   if(h->ecat_calibration_fctr>0.0) h->quant_scale*=h->ecat_calibration_fctr;
00546   fptr=_fdata; mptr=mdata;
00547   if(h->data_type==BYTE_TYPE) {
00548     for(i=0; i<pxlNr; i++, mptr++, fptr++)
00549       *fptr=h->quant_scale*(float)(*mptr);
00550   } else if(h->data_type==VAX_I2 || h->data_type==SUN_I2) {
00551     for(i=0; i<pxlNr; i++, mptr+=2, fptr++) {
00552       sptr=(short int*)mptr;
00553       *fptr=h->quant_scale*(float)(*sptr);
00554     }
00555   } else if(h->data_type==VAX_I4 || h->data_type==SUN_I4) {
00556     for(i=0; i<pxlNr; i++, mptr+=4, fptr++) {
00557       iptr=(int*)mptr;
00558       *fptr=h->quant_scale*(float)(*iptr);
00559     }
00560   } else if(h->data_type==VAX_R4 || h->data_type==IEEE_R4) {
00561     memcpy(fptr, mptr, pxlNr*4);
00562     for(i=0; i<pxlNr; i++, fptr++) *fptr *= h->quant_scale;
00563   }
00564   free(mdata);
00565   *fdata=_fdata;
00566 
00567   return(0);
00568 }
00569 /*****************************************************************************/
00570 
00571 /*****************************************************************************/
00586 int ecat63ReadScanMatrix(FILE *fp, int first_block, int last_block, ECAT63_scanheader *h, float **fdata) {
00587   int i, ret, blockNr, pxlNr;
00588   char *mdata, *mptr;
00589   float *_fdata, *fptr;
00590   short int *sptr;
00591   int *iptr;
00592   
00593   
00594   if(ECAT63_TEST) printf("ecat63ReadScanMatrix(fp, %d, %d, hdr, fdata)\n",
00595     first_block, last_block);
00596   if(fp==NULL || first_block<=MatFirstDirBlk || h==NULL) {
00597     sprintf(ecat63errmsg, "invalid function parameter.\n");
00598     return(1);
00599   }
00600   *fdata=(float*)NULL;
00601   
00602   /* Read subheader */
00603   ret=ecat63ReadScanheader(fp, first_block, h);
00604   if(ret) {
00605     sprintf(ecat63errmsg, "cannot read subheader (%d).\n", ret);
00606     return(5);
00607   }
00608   if(ECAT63_TEST>4) ecat63PrintScanheader(h, stdout);
00609   pxlNr=h->dimension_1*h->dimension_2;
00610   if(pxlNr<=0) {
00611     sprintf(ecat63errmsg, "invalid matrix dimension.\n");
00612     return(6);  
00613   }
00614   
00615   /* Read matrix data */
00616   blockNr=last_block-first_block; if(blockNr<1) return(0);
00617   mdata=(char*)malloc(blockNr*MatBLKSIZE);
00618   if(mdata==NULL) {
00619     sprintf(ecat63errmsg, "cannot allocate memory.\n");
00620     return(8);  
00621   }
00622   mptr=mdata;
00623   ret=ecat63ReadMatdata(fp, first_block+1, blockNr, mptr, h->data_type);
00624   if(ret || mdata==NULL) {
00625     sprintf(ecat63errmsg, "cannot read matrix data (%d).\n", ret);
00626     free(mdata); return(9);
00627   }
00628   
00629   /* Allocate memory for float data */
00630   _fdata=(float*)malloc(pxlNr*sizeof(float));
00631   if(_fdata==NULL) {
00632     sprintf(ecat63errmsg, "cannot allocate memory.\n");
00633     free(mdata); return(11);  
00634   }
00635 
00636   /* Convert matrix data to floats */
00637   fptr=_fdata; mptr=mdata;
00638   if(h->data_type==BYTE_TYPE) {
00639     for(i=0; i<pxlNr; i++, mptr++, fptr++)
00640       *fptr=h->scale_factor*(float)(*mptr);
00641   } else if(h->data_type==VAX_I2 || h->data_type==SUN_I2) {
00642     for(i=0; i<pxlNr; i++, mptr+=2, fptr++) {
00643       sptr=(short int*)mptr;
00644       *fptr=h->scale_factor*(float)(*sptr);
00645     }
00646   } else if(h->data_type==VAX_I4 || h->data_type==SUN_I4) {
00647     for(i=0; i<pxlNr; i++, mptr+=4, fptr++) {
00648       iptr=(int*)mptr;
00649       *fptr=h->scale_factor*(float)(*iptr);
00650     }
00651   } else if(h->data_type==VAX_R4 || h->data_type==IEEE_R4) {
00652     memcpy(fptr, mptr, pxlNr*4);
00653     for(i=0; i<pxlNr; i++, fptr++) *fptr *= h->scale_factor;
00654   }
00655   free(mdata);
00656   *fdata=_fdata;
00657 
00658   return(0);
00659 }
00660 /*****************************************************************************/
00661 
00662 /*****************************************************************************/
00671 float ecat63rFloat(void *bufi, int isvax, int islittle) {
00672   union {unsigned int ul; float f;} t;
00673 
00674   memcpy(&t.ul, bufi, 4); if(t.ul==0) {return(0.0);}
00675   if(isvax) { /* if input is in VAX format */
00676     /* Swap words on i386 and bytes on SUN */
00677     if(islittle) swawip(&t.ul, 4); else swabip(&t.ul, 4);
00678     t.ul-=(2L<<23); /* subtract 2 from exp */
00679   } else { /* input is in i386 format */
00680     if(!islittle) swawbip(&t.ul, 4); /* Switch words and bytes on SUN */
00681   }
00682   return(t.f);
00683 }
00684 
00694 int ecat63rInt(void *bufi, int isvax, int islittle) {
00695   int i;
00696 
00697   /* Swap both words and bytes on SUN */
00698   memcpy(&i, bufi, 4); if(!islittle) swawbip(&i, 4);
00699   return(i);
00700 }
00701 /*****************************************************************************/
00702 
00703 /*****************************************************************************/
00711 int ecat63pxlbytes(short int data_type) {
00712   int byteNr=0;
00713   switch(data_type) {
00714     case BYTE_TYPE: byteNr=1; break;
00715     case VAX_I2:
00716     case SUN_I2: byteNr=2; break;
00717     case VAX_I4:
00718     case VAX_R4:
00719     case IEEE_R4:
00720     case SUN_I4: byteNr=4; break;
00721   }
00722   return(byteNr);
00723 }
00724 /*****************************************************************************/
00725 
00726 /*****************************************************************************/
00727