ecat63p.c
Go to the documentation of this file.
00001 /******************************************************************************
00002 
00003   Copyright (c) 2003-2008 Turku PET Centre
00004 
00005   Library:     ecat63p.c
00006   Description: Procedures for printing ECAT 6.3 header contents.
00007 
00008   This library is free software; you can redistribute it and/or
00009   modify it under the terms of the GNU Lesser General Public
00010   License as published by the Free Software Foundation; either
00011   version 2.1 of the License, or (at your option) any later version.
00012 
00013   This library is distributed in the hope that it will be useful,
00014   but WITHOUT ANY WARRANTY; without even the implied warranty of
00015   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
00016   See the GNU Lesser General Public License for more details:
00017   http://www.gnu.org/copyleft/lesser.html
00018 
00019   You should have received a copy of the GNU Lesser General Public License
00020   along with this library/program; if not, write to the Free Software
00021   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA 
00022 
00023   Turku PET Centre, Turku, Finland, http://www.turkupetcentre.fi/
00024 
00025   Modification history:
00026   2003-07-21 Vesa Oikonen
00027     Contents separated from ecat63.c.
00028   2004-02-08 VO
00029     ecat63PrintMainheader() prints also sw_version.
00030   2004-09-20 VO
00031     Doxygen style comments.
00032   2007-03-12 VO
00033     Header print functions changed: file pointer is required as second
00034     argument, set to stdout to make it work as before;
00035     header contents are now printed one value per line.
00036   2008-07-24 VO
00037     Added ecat63PrintSubheader().
00038 
00039 
00040 ******************************************************************************/
00041 #include <stdio.h>
00042 #include <stdlib.h>
00043 #include <math.h>
00044 #include <ctype.h>
00045 #include <string.h>
00046 #include <unistd.h>
00047 #include <time.h>
00048 /*****************************************************************************/
00049 #include "swap.h"
00050 #include "include/ecat63.h"
00051 /*****************************************************************************/
00052 
00053 /*****************************************************************************/
00054 
00055 /*****************************************************************************/
00062 void ecat63PrintMainheader(ECAT63_mainheader *h, FILE *fp) {
00063   char tmp[32];
00064 
00065   if(ECAT63_TEST) fprintf(stdout, "ecat63PrintMainheader()\n");
00066   fprintf(fp, "original_file_name := %.20s\n", h->original_file_name);
00067   fprintf(fp, "data_type := %d (%s)\n", h->data_type, ecat63Datatype(h->data_type));
00068   fprintf(fp, "system_type := %d\n", h->system_type);
00069   if(h->file_type==1) strcpy(tmp, "sinogram");
00070   else if(h->file_type==2) strcpy(tmp, "image");
00071   else if(h->file_type==3) strcpy(tmp, "attenuation");
00072   else if(h->file_type==4) strcpy(tmp, "normalization");
00073   else strcpy(tmp, "unknown");
00074   fprintf(fp, "file_type := %d (%s)\n", h->file_type, tmp);
00075   fprintf(fp, "sw_version := %d\n", h->sw_version);
00076   fprintf(fp, "scan start := %02d.%02d.%04d %02d:%02d:%02d\n",
00077     h->scan_start_day, h->scan_start_month, h->scan_start_year,
00078     h->scan_start_hour, h->scan_start_minute, h->scan_start_second);
00079   fprintf(fp, "isotope_code := %.8s\n", h->isotope_code);
00080   fprintf(fp, "isotope_halflife := %E sec\n", h->isotope_halflife);
00081   fprintf(fp, "radiopharmaceutical := %.32s\n", h->radiopharmaceutical);
00082   fprintf(fp, "gantry_tilt := %g\n", h->gantry_tilt);
00083   fprintf(fp, "gantry_rotation := %g\n", h->gantry_rotation);
00084   fprintf(fp, "bed_elevation := %g\n", h->bed_elevation);
00085   fprintf(fp, "axial_fov := %g\n", h->axial_fov);
00086   fprintf(fp, "transaxial_fov := %g\n", h->transaxial_fov);
00087   fprintf(fp, "calibration_factor := %g\n", h->calibration_factor);
00088   fprintf(fp, "calibration_units := %d (%s)\n", h->calibration_units, ecat63Unit(h->calibration_units));
00089   fprintf(fp, "study_name := %.12s\n", h->study_name);
00090   fprintf(fp, "patient_id := %.32s\n", h->patient_id);
00091   fprintf(fp, "patient_name := %.32s\n", h->patient_name);
00092   fprintf(fp, "sex := %c\nage := %.10s\nheight := %.10s\nweigth := %.10s\ndexterity := %c\n",
00093     h->patient_sex, h->patient_age, h->patient_height, h->patient_weight,
00094     h->patient_dexterity);
00095   fprintf(fp, "physician_name := %.32s\noperator_name := %.32s\n",
00096     h->physician_name, h->operator_name);
00097   fprintf(fp, "study_description := %.32s\n", h->study_description);
00098   fprintf(fp, "num_planes := %d\n", h->num_planes);
00099   fprintf(fp, "num_frames := %d\n", h->num_frames);
00100   fprintf(fp, "num_gates := %d\n", h->num_gates);
00101   fprintf(fp, "num_bed_pos := %d\n", h->num_bed_pos);
00102   fprintf(fp, "init_bed_position := %g\n", h->init_bed_position);
00103   fprintf(fp, "plane_separation := %g cm\n", h->plane_separation);
00104   fprintf(fp, "user_process_code := %.10s\n", h->user_process_code);
00105 }
00106 /*****************************************************************************/
00107 
00108 /*****************************************************************************/
00115 void ecat63PrintImageheader(ECAT63_imageheader *h, FILE *fp) {
00116   int i;
00117 
00118   if(ECAT63_TEST) printf("ecat63PrintImageheader()\n");
00119   fprintf(fp, "data_type := %d (%s)\n", h->data_type, ecat63Datatype(h->data_type));
00120   fprintf(fp, "dimension_1 := %d\n", h->dimension_1);
00121   fprintf(fp, "dimension_2 := %d\n", h->dimension_2);
00122   fprintf(fp, "x_origin := %g\ny_origin := %g\nrecon_scale := %g\n",
00123     h->x_origin, h->y_origin, h->recon_scale);
00124   fprintf(fp, "quant_scale := %g\nimage_min := %d\nimage_max := %d\n",
00125     h->quant_scale, h->image_min, h->image_max);
00126   fprintf(fp, "slice_width := %g\npixel_size := %g\n", h->slice_width, h->pixel_size);
00127   fprintf(fp, "frame_start_time := %d\nframe_duration := %d\n",
00128     h->frame_start_time, h->frame_duration );
00129   fprintf(fp, "reconstruction_start := %02d.%02d.%04d %02d:%02d:%02d\n",
00130     h->recon_start_day, h->recon_start_month, h->recon_start_year,
00131     h->recon_start_hour, h->recon_start_min, h->recon_start_sec);
00132   fprintf(fp, "filter_code := %d\nimage_rotation := %g\nintrinsic_tilt := %g\n",
00133     h->filter_code, h->image_rotation, h->intrinsic_tilt);
00134   fprintf(fp, "filter_params :=");
00135   for(i=0; i<6; i++) fprintf(fp, " %g", h->filter_params[i]); fprintf(fp, "\n");
00136   fprintf(fp, "plane_eff_corr_fctr := %g\ndecay_corr_fctr := %g\nloss_corr_fctr := %g\n",
00137     h->plane_eff_corr_fctr, h->decay_corr_fctr, h->loss_corr_fctr);
00138   fprintf(fp, "quant_units := %d (%s)\n", h->quant_units, ecat63Unit(h->quant_units));
00139   fprintf(fp, "ecat_calibration_fctr := %g\nwell_counter_cal_fctr := %g\n",
00140     h->ecat_calibration_fctr, h->well_counter_cal_fctr);
00141   fprintf(fp, "annotation := %.40s\n", h->annotation);
00142 }
00143 /*****************************************************************************/
00144 
00145 /*****************************************************************************/
00152 void ecat63PrintScanheader(ECAT63_scanheader *h, FILE *fp) {
00153   int i;
00154 
00155   if(ECAT63_TEST) printf("ecat63PrintScanheader()\n");
00156   fprintf(fp, "data_type := %d (%s)\n", h->data_type, ecat63Datatype(h->data_type));
00157   fprintf(fp, "dimension_1 := %d\n", h->dimension_1);
00158   fprintf(fp, "dimension_2 := %d\n", h->dimension_2);
00159   fprintf(fp, "sample_distance := %g cm\n", h->sample_distance);
00160   fprintf(fp, "isotope_halflife := %g sec\n", h->isotope_halflife);
00161   fprintf(fp, "gate_duration := %d\nr_wave_offset := %d\n",
00162     h->gate_duration, h->r_wave_offset);
00163   fprintf(fp, "scale_factor := %g\n", h->scale_factor);
00164   fprintf(fp, "scan_min := %d\nscan_max := %d\n", h->scan_min, h->scan_max);
00165   fprintf(fp, "prompts := %d\ndelayed := %d\nmultiples := %d\nnet_trues := %d\n",
00166     h->prompts, h->delayed, h->multiples, h->net_trues);
00167   fprintf(fp, "cor_singles :=");
00168   for(i=0; i<16; i++) fprintf(fp, " %8.0f", h->cor_singles[i]); printf("\n");
00169   fprintf(fp, "uncor_singles :=");
00170   for(i=0; i<16; i++) fprintf(fp, " %8.0f", h->uncor_singles[i]); printf("\n");
00171   fprintf(fp, "tot_avg_cor := %g\ntot_avg_uncor := %g\n", h->tot_avg_cor, h->tot_avg_uncor);
00172   fprintf(fp, "total_coin_rate := %d\n", h->total_coin_rate);
00173   fprintf(fp, "frame_start_time := %d\nframe_duration := %d\n",
00174     h->frame_start_time, h->frame_duration);
00175   fprintf(fp, "loss_correction_fctr := %g\n", h->loss_correction_fctr);
00176 }
00177 /*****************************************************************************/
00178 
00179 /*****************************************************************************/
00186 void ecat63PrintAttnheader(ECAT63_attnheader *h,  FILE *fp) {
00187   if(ECAT63_TEST) printf("ecat63PrintAttnheader()\n");
00188   fprintf(fp, "data_type := %d (%s)\n", h->data_type, ecat63Datatype(h->data_type));
00189   fprintf(fp, "dimension_1 := %d\n", h->dimension_1);
00190   fprintf(fp, "dimension_2 := %d\n", h->dimension_2);
00191   fprintf(fp, "sample_distance := %g cm\n", h->sample_distance);
00192   fprintf(fp, "attenuation_type := %d\n", h->attenuation_type);
00193   fprintf(fp, "scale_factor := %g\n", h->scale_factor);
00194   fprintf(fp, "x_origin := %g\ny_origin := %g\nx_radius := %g\ny_radius := %g\n",
00195     h->x_origin, h->y_origin, h->x_radius, h->y_radius);
00196   fprintf(fp, "tilt_angle := %g\nattenuation_coeff := %g\n",
00197     h->tilt_angle, h->attenuation_coeff);
00198 }
00199 /*****************************************************************************/
00200 
00201 /*****************************************************************************/
00208 void ecat63PrintNormheader(ECAT63_normheader *h, FILE *fp) {
00209   if(ECAT63_TEST) printf("ecat63PrintNormheader()\n");
00210   fprintf(fp, "data_type := %d (%s)\n", h->data_type, ecat63Datatype(h->data_type));
00211   fprintf(fp, "dimension_1 := %d\n", h->dimension_1);
00212   fprintf(fp, "dimension_2 := %d\n", h->dimension_2);
00213   fprintf(fp, "scale_factor := %g\n", h->scale_factor);
00214   fprintf(fp, "norm time := %02d.%02d.%04d %02d:%02d:%02d\n",
00215     h->norm_day, h->norm_month, h->norm_year,
00216     h->norm_hour, h->norm_minute, h->norm_second);
00217 }
00218 /*****************************************************************************/
00219 
00220 /*****************************************************************************/
00227 char *ecat63Datatype(short int dtype) {
00228   static char *ecat63_datatype[]={
00229   /*  0 */  "Unknown",
00230   /*  1 */  "BYTE_TYPE",
00231   /*  2 */  "VAX_I2",
00232   /*  3 */  "VAX_I4",
00233   /*  4 */  "VAX_R4",
00234   /*  5 */  "IEEE_R4",
00235   /*  6 */  "SUN_I2",
00236   /*  7 */  "SUN_I4",
00237   /*  8 */  "Unknown",
00238   /*  9 */  "Unknown"
00239   };
00240   if(dtype>=0 && dtype<10) return(ecat63_datatype[dtype]);
00241   else return(ecat63_datatype[0]);
00242 }
00243 
00250 char *ecat63Unit(short int dunit) {
00251   static char *ecat63_unit[]={
00252   /*  0 */  "Unknown",
00253   /*  1 */  "Unknown",
00254   /*  2 */  "ECAT counts",
00255   /*  3 */  "uCi/ml",
00256   /*  4 */  "LMRGlu",
00257   /*  5 */  "LMRUGlu umol/min/100g",
00258   /*  6 */  "LMRUGlu mg/min/100g",
00259   /*  7 */  "nCi/mL",
00260   /*  8 */  "Well counts",
00261   /*  9 */  "Becquerels",
00262   /* 10 */  "kBq/mL",
00263   /* 11 */  "1/min",
00264   /* 12 */  "mL/min/100g",
00265   /* 13 */  "sec*kBq/mL",
00266   /* 14 */  "sec*nCi/mL",
00267   /* 15 */  "1/sec",
00268   /* 16 */  "Unitless",
00269   /* 17 */  "Unknown"
00270   };
00271   if(dunit>=0 && dunit<18) return(ecat63_unit[dunit]);
00272   else return(ecat63_unit[0]);
00273 }
00274 /*****************************************************************************/
00282 void float2parts(float *buf) {
00283   unsigned int u, exp, mantissa;
00284   char sign;
00285 
00286   memcpy(&u, buf, 4); if(u & 1L<<31) sign='-'; else sign='+';
00287   exp=u<<1; exp=exp>>24; mantissa=u<<9; mantissa=mantissa>>9;
00288   printf("%e = %c (%u/8388608 + 1)*2^(%u-127)\n", *buf, sign, mantissa, exp);
00289 }
00290 /*****************************************************************************/
00291 
00292 /*****************************************************************************/
00296 int ecat6PrintSubheader(
00298   ECAT63_mainheader mh,
00300   FILE *fp,
00302   int plane,
00304   int frame,
00306   FILE *ofp
00307 ) {
00308   int                 mi, ret, nr=0;
00309   static MATRIXLIST   mlist;
00310   ECAT63_imageheader  image_header;
00311   ECAT63_scanheader   scan_header;
00312   ECAT63_attnheader   attn_header;
00313   ECAT63_normheader   norm_header;
00314   Matval              matval;
00315 
00316 
00317   /*  Read matrix list and nr */
00318   ecat63InitMatlist(&mlist);
00319   ret=ecat63ReadMatlist(fp, &mlist);
00320   if(ret) {
00321     fprintf(stderr, "Error (%d): cannot read matrix list.\n", ret);
00322     return 2;
00323   }
00324   if(mlist.matrixNr<=0) {
00325     fprintf(stderr, "Error: matrix list is empty.\n");
00326     return 2;
00327   }
00328   if(ECAT63_TEST>1) ecat63PrintMatlist(&mlist);
00329 
00330   /*
00331    *  Read and print subheaders one at a time
00332    */
00333   for(mi=nr=0; mi<mlist.matrixNr; mi++) {
00334     /* Get plane and frame nr */
00335     mat_numdoc(mlist.matdir[mi].matnum, &matval);
00336     /* Check if this is supposed to be listed or not */
00337     if(frame>=0 && frame!=matval.frame) continue;
00338     if(plane>=0 && plane!=matval.plane) continue;
00339     /* Read subheader */
00340     if(mh.file_type==IMAGE_DATA)
00341       ret=ecat63ReadImageheader(fp, mlist.matdir[mi].strtblk, &image_header);
00342     else if(mh.file_type==RAW_DATA)
00343       ret=ecat63ReadScanheader(fp, mlist.matdir[mi].strtblk, &scan_header);
00344     else if(mh.file_type==ATTN_DATA)
00345       ret=ecat63ReadAttnheader(fp, mlist.matdir[mi].strtblk, &attn_header);
00346     else if(mh.file_type==NORM_DATA)
00347       ret=ecat63ReadNormheader(fp, mlist.matdir[mi].strtblk, &norm_header);
00348     if(ret) {
00349       fprintf(stderr, "Error: cannot read matrix %u subheader.\n",
00350          mlist.matdir[mi].matnum);
00351       ecat63EmptyMatlist(&mlist); return 4;
00352     }
00353     /* Print subheader */
00354     fprintf(fp, "Matrix: plane %d frame %d gate %d bed %d\n",
00355       matval.plane, matval.frame, matval.gate, matval.bed);
00356     if(mh.file_type==IMAGE_DATA)
00357       ecat63PrintImageheader(&image_header, ofp);
00358     else if(mh.file_type==RAW_DATA)
00359       ecat63PrintScanheader(&scan_header, ofp);
00360     else if(mh.file_type==ATTN_DATA)
00361       ecat63PrintAttnheader(&attn_header, ofp);
00362     else if(mh.file_type==NORM_DATA)
00363       ecat63PrintNormheader(&norm_header, ofp);
00364     nr++; // counter
00365   } /* next matrix */
00366   ecat63EmptyMatlist(&mlist);
00367   
00368   if(nr==0 && (plane>=0 || frame>=0)) {
00369     fprintf(stderr, "Error: specified matrices not found.\n");
00370     return(11);
00371   }
00372 
00373   return(0);
00374 }
00375 /*****************************************************************************/
00376 
00377 /*****************************************************************************/
00378