TPCCLIB
Loading...
Searching...
No Matches
img_e63.c File Reference

ECAT 6.3 I/O routines for IMG data. More...

#include "libtpcimgio.h"

Go to the source code of this file.

Functions

int ecat63ReadAllToImg (const char *fname, IMG *img)
int ecat63WriteAllImg (const char *fname, IMG *img)
int ecat63ReadPlaneToImg (const char *fname, IMG *img)
int ecat63AddImg (const char *fname, IMG *img)
int imgEcat63Supported (ECAT63_mainheader *h)
void imgGetEcat63MHeader (IMG *img, ECAT63_mainheader *h)
void imgSetEcat63MHeader (IMG *img, ECAT63_mainheader *h)
int imgGetEcat63Fileformat (ECAT63_mainheader *h)
int imgReadEcat63Header (const char *fname, IMG *img)
int imgReadEcat63FirstFrame (const char *fname, IMG *img)
int imgReadEcat63Frame (const char *fname, int frame_to_read, IMG *img, int frame_index)
int imgWriteEcat63Frame (const char *fname, int frame_to_write, IMG *img, int frame_index)
void imgSetEcat63SHeader (IMG *img, void *h)

Detailed Description

ECAT 6.3 I/O routines for IMG data.

Author
Vesa Oikonen

Definition in file img_e63.c.

Function Documentation

◆ ecat63AddImg()

int ecat63AddImg ( const char * fname,
IMG * img )

Adds all matrices in memory to the ECAT file. If ECAT file does not exist, it is created.

Please note that existing ECAT file is NOT saved as bak file.

Returns
0 if ok, 1 invalid input, 2 image status is not 'occupied', 3 failed to open file for reading, 4 failed to allocate memory for data, 9 failed to write data, 21 invalid matrix list, 22 failed to write main header
Parameters
fnameName of the output ECAT 6.3 file.
imgPointer to data structure from which the data is written.

Definition at line 878 of file img_e63.c.

883 {
884 int n, m, ret=0, add=0;
885 int frameNr, planeNr;
886 int frame, plane, prev_plane;
887 float f, fmax, fmin, g, scale;
888 short int *sdata, *sptr, smin, smax;
889 FILE *fp;
890 ECAT63_mainheader main_header;
891 ECAT63_imageheader image_header;
892 ECAT63_scanheader scan_header;
893 MATRIXLIST mlist;
894 MatDir tmpMatdir;
895 Matval matval;
896 struct tm tm;
897
898 if(IMG_TEST) printf("ecat63AddImg(%s, *img)\n", fname);
899 if(IMG_TEST>4) imgInfo(img);
900 /* Check the arguments */
901 if(fname==NULL) {strcpy(ecat63errmsg, "invalid ECAT filename"); return(1);}
902 if(img==NULL || img->status!=IMG_STATUS_OCCUPIED) {
903 strcpy(ecat63errmsg, "image data is empty"); return(2);}
904 if(img->_dataType<1) img->_dataType=VAX_I2;
905
906 /* Check image size */
907 if(img->dimt>511 || img->dimz>255 || img->dimx>1024 || img->dimy>1024) {
908 strcpy(ecat63errmsg, "too large matrix size"); return(2);}
909
910 /* Initiate headers */
911 memset(&main_header, 0, sizeof(ECAT63_mainheader));
912 memset(&image_header,0, sizeof(ECAT63_imageheader));
913 memset(&scan_header, 0, sizeof(ECAT63_scanheader));
914
915 /* Make a main header */
916 main_header.sw_version=2;
917 main_header.num_planes=img->dimz;
918 main_header.num_frames=img->dimt;
919 main_header.num_gates=1;
920 main_header.num_bed_pos=1;
921 if(img->type==IMG_TYPE_IMAGE) main_header.file_type=IMAGE_DATA;
922 else if(img->type==IMG_TYPE_RAW) main_header.file_type=RAW_DATA;
923 else {strcpy(ecat63errmsg, "invalid filetype"); return(1);}
924 main_header.data_type=img->_dataType;
925 if(img->scanner>0) main_header.system_type=img->scanner;
927 main_header.calibration_factor=1.0;
928 main_header.axial_fov=img->axialFOV/10.0;
929 main_header.transaxial_fov=img->transaxialFOV/10.0;
930 main_header.plane_separation=img->sizez/10.0;
931 main_header.calibration_units=imgUnitToEcat6(img);
932 strncpy(main_header.radiopharmaceutical, img->radiopharmaceutical, 32);
933 if(gmtime_r((time_t*)&img->scanStart, &tm)!=NULL) {
934 main_header.scan_start_year=tm.tm_year+1900;
935 main_header.scan_start_month=tm.tm_mon+1;
936 main_header.scan_start_day=tm.tm_mday;
937 main_header.scan_start_hour=tm.tm_hour;
938 main_header.scan_start_minute=tm.tm_min;
939 main_header.scan_start_second=tm.tm_sec;
940 if(IMG_TEST>2) {
941 printf(" img->scanStart := %ld\n", (long int)img->scanStart);
942 printf(" -> tm_year := %d\n", tm.tm_year);
943 printf(" -> tm_hour := %d\n", tm.tm_hour);
944 }
945 } else {
946 main_header.scan_start_year=1900;
947 main_header.scan_start_month=1;
948 main_header.scan_start_day=1;
949 main_header.scan_start_hour=0;
950 main_header.scan_start_minute=0;
951 main_header.scan_start_second=0;
952 if(IMG_TEST>0) printf("invalid scan_start_time in IMG\n");
953 }
954 main_header.isotope_halflife=img->isotopeHalflife;
955 strcpy(main_header.isotope_code, imgIsotope(img));
956 strlcpy(main_header.study_name, img->studyNr, 12);
957 strcpy(main_header.patient_name, img->patientName);
958 strcpy(main_header.patient_id, img->patientID);
959 strlcpy(main_header.study_description, img->studyDescription, 32);
960 strlcpy(main_header.user_process_code, img->userProcessCode, 10);
961
962 /* Check if the ECAT file exists */
963 if(access(fname, 0) != -1) {
964 if(IMG_TEST) printf("Opening existing ECAT file.\n");
965 add=1;
966 /* Open the ECAT file */
967 if((fp=fopen(fname, "rb+")) == NULL) {
968 strcpy(ecat63errmsg, "cannot create ECAT file"); return(3);}
969 /* Read main header */
970 if((ret=ecat63ReadMainheader(fp, &main_header))) {
971 sprintf(ecat63errmsg, "cannot read main header (%d)", ret);
972 fclose(fp); return(3);
973 }
974 fflush(fp);
975 /* Check that filetypes are matching */
976 if((main_header.file_type==IMAGE_DATA && img->type==IMG_TYPE_IMAGE) ||
977 (main_header.file_type==RAW_DATA && img->type==IMG_TYPE_RAW)) {
978 /* ok */
979 } else {
980 sprintf(ecat63errmsg, "cannot add different filetype");
981 fclose(fp); return(3);
982 }
983 } else {
984 if(IMG_TEST) printf("ECAT file does not exist.\n");
985 add=0;
986 /* Create output ECAT file */
987 fp=ecat63Create(fname, &main_header);
988 if(fp==NULL) {strcpy(ecat63errmsg, "cannot create ECAT file"); return(3);}
989 }
990 if(IMG_TEST) ecat63PrintMainheader(&main_header, stdout);
991
992 /* Allocate memory for matrix data */
993 sdata=(short int*)malloc((size_t)img->dimx*img->dimy*sizeof(short int));
994 if(sdata==NULL) {strcpy(ecat63errmsg, "out of memory"); return(4);}
995
996 /* Set the subheader contents, as far as possible */
997 switch(main_header.file_type) {
998 case RAW_DATA:
999 scan_header.data_type=main_header.data_type;
1000 scan_header.dimension_1=img->dimx;
1001 scan_header.dimension_2=img->dimy;
1002 scan_header.frame_duration_sec=1;
1003 scan_header.scale_factor=1.0;
1004 scan_header.frame_start_time=0;
1005 scan_header.frame_duration=1000;
1006 scan_header.loss_correction_fctr=1.0;
1007 scan_header.sample_distance=(img->sampleDistance)/10.0;
1008 /*if(IMG_TEST) ecat63PrintScanheader(&scan_header);*/
1009 break;
1010 case IMAGE_DATA:
1011 image_header.data_type=main_header.data_type;
1012 image_header.num_dimensions=2;
1013 image_header.dimension_1=img->dimx;
1014 image_header.dimension_2=img->dimy;
1015 image_header.recon_scale=img->zoom;
1016 image_header.quant_scale=1.0;
1017 image_header.slice_width=img->sizez/10.;
1018 image_header.pixel_size=img->sizex/10.;
1019 image_header.frame_start_time=0;
1020 image_header.frame_duration=1000;
1021 image_header.plane_eff_corr_fctr=1.0;
1022 image_header.decay_corr_fctr=0.0;
1023 image_header.loss_corr_fctr=1.0;
1024 image_header.ecat_calibration_fctr=1.0;
1025 image_header.well_counter_cal_fctr=1.0;
1026 image_header.quant_units=main_header.calibration_units;
1027 /*if(IMG_TEST) ecat63PrintImageheader(&image_header);*/
1028 break;
1029 }
1030
1031 /* Write one matrix at a time */
1032 n=0;
1033 for(plane=1; plane<=img->dimz;plane++) for(frame=1;frame<=img->dimt;frame++) {
1034 /* Scale data to short ints */
1035 /* Search min and max */
1036 fmin=fmax=f=img->m[plane-1][0][0][frame-1];
1037 for(int i=0; i<img->dimy; i++) for(int j=0; j<img->dimx; j++) {
1038 f=img->m[plane-1][i][j][frame-1];
1039 if(f>fmax) fmax=f; else if(f<fmin) fmin=f;
1040 }
1041 /* Calculate scaling factor */
1042 if(fabs(fmin)>fabs(fmax)) g=fabs(fmin); else g=fabs(fmax);
1043 if(!isfinite(g)) {
1044 sprintf(ecat63errmsg, "invalid pixel values on pl%02d fr%02d (%d).", plane, frame, ret);
1045 fclose(fp); remove(fname); free(sdata);
1046 return(1);
1047 }
1048 if(g>0.0) scale=32766./g; else scale=1.0;
1049 if(!isfinite(scale)) scale=1.0;
1050 /*printf("fmin=%e fmax=%e g=%e scale=%e\n", fmin, fmax, g, scale);*/
1051 /* Scale matrix data to shorts */
1052 sptr=sdata;
1053 for(int i=0; i<img->dimy; i++) for(int j=0; j<img->dimx; j++) {
1054 *sptr=(short int)temp_roundf(scale*img->m[plane-1][i][j][frame-1]);
1055 sptr++;
1056 }
1057 /* Calculate and set subheader min&max and scale */
1058 smin=(short int)temp_roundf(scale*fmin);
1059 smax=(short int)temp_roundf(scale*fmax);
1060 if(main_header.file_type==RAW_DATA) {
1061 scan_header.scan_min=smin; scan_header.scan_max=smax;
1062 scan_header.scale_factor=1.0/scale;
1063 scan_header.frame_start_time=(int)temp_roundf(1000.*img->start[frame-1]);
1064 scan_header.frame_duration=
1065 (int)temp_roundf(1000.*(img->end[frame-1]-img->start[frame-1]));
1066 scan_header.prompts=temp_roundf(img->prompts[frame-1]);
1067 scan_header.delayed=temp_roundf(img->randoms[frame-1]);
1068 } else if(main_header.file_type==IMAGE_DATA) {
1069 image_header.image_min=smin; image_header.image_max=smax;
1070 image_header.quant_scale=1.0/scale;
1071 image_header.frame_start_time=(int)temp_roundf(1000.*img->start[frame-1]);
1072 image_header.frame_duration=
1073 (int)temp_roundf(1000.*(img->end[frame-1]-img->start[frame-1]));
1074 /* Set decay correction factor */
1076 image_header.decay_corr_fctr=img->decayCorrFactor[frame-1];
1077 else
1078 image_header.decay_corr_fctr=0.0;
1079 }
1080 /* Write matrix data */
1081 m=mat_numcod(frame, img->planeNumber[plane-1], 1, 0, 0);
1082 sptr=sdata;
1083 if(IMG_TEST) printf(" writing matnum=%d\n", m);
1084 if(main_header.file_type==RAW_DATA) {
1085 /*if(IMG_TEST) ecat63PrintScanheader(&scan_header);*/
1086 ret=ecat63WriteScan(fp, m, &scan_header, sptr);
1087 } else if(main_header.file_type==IMAGE_DATA) {
1088 /*if(IMG_TEST) ecat63PrintImageheader(&image_header);*/
1089 ret=ecat63WriteImage(fp, m, &image_header, sptr);
1090 }
1091 if(ret) {
1092 sprintf(ecat63errmsg, "cannot write data on pl%02d fr%02d (%d).",
1093 plane, frame, ret);
1094 fclose(fp); remove(fname); free(sdata);
1095 return(9);
1096 }
1097 n++;
1098 } /* next matrix */
1099 if(IMG_TEST) printf(" %d matrices written in %s\n", n, fname);
1100
1101 /* Free matrix memory */
1102 free(sdata);
1103
1104 /* If matrices were added, set main header */
1105 if(add==1) {
1106 fflush(fp);
1107 /* Read matrix list */
1108 ecat63InitMatlist(&mlist);
1109 ret=ecat63ReadMatlist(fp, &mlist, IMG_TEST);
1110 if(ret) {
1111 sprintf(ecat63errmsg, "cannot read matrix list (%d)", ret);
1112 fclose(fp); return(21);
1113 }
1114 if(mlist.matrixNr<=0) {
1115 strcpy(ecat63errmsg, "matrix list is empty"); fclose(fp); return(21);}
1116 /* Sort matrix list */
1117 for(int i=0; i<mlist.matrixNr-1; i++) for(int j=i+1; j<mlist.matrixNr; j++) {
1118 if(mlist.matdir[i].matnum>mlist.matdir[j].matnum) {
1119 tmpMatdir=mlist.matdir[i];
1120 mlist.matdir[i]=mlist.matdir[j]; mlist.matdir[j]=tmpMatdir;
1121 }
1122 }
1123 /* Calculate the number of planes and frames */
1124 prev_plane=plane=-1; /*prev_frame=frame=-1;*/ frameNr=planeNr=0;
1125 for(int m=0; m<mlist.matrixNr; m++) {
1126 mat_numdoc(mlist.matdir[m].matnum, &matval);
1127 plane=matval.plane; frame=matval.frame;
1128 if(plane!=prev_plane) {frameNr=1; planeNr++;} else {frameNr++;}
1129 prev_plane=plane; /*prev_frame=frame;*/
1130 } /* next matrix */
1131 ecat63EmptyMatlist(&mlist);
1132 main_header.num_planes=planeNr; main_header.num_frames=frameNr;
1133 /* Write main header */
1134 ret=ecat63WriteMainheader(fp, &main_header);
1135 if(ret) {
1136 sprintf(ecat63errmsg, "cannot write mainheader (%d)", ret);
1137 fclose(fp); return(22);
1138 }
1139 }
1140
1141 /* Close file and free memory */
1142 fclose(fp);
1143
1144 return(0);
1145}
struct tm * gmtime_r(const time_t *t, struct tm *tm)
Convert time_t to GMT struct tm.
Definition datetime.c:22
char ecat63errmsg[128]
Definition ecat63h.c:7
int ecat63ReadMatlist(FILE *fp, MATRIXLIST *ml, int verbose)
Definition ecat63ml.c:46
void ecat63InitMatlist(MATRIXLIST *mlist)
Definition ecat63ml.c:20
void ecat63EmptyMatlist(MATRIXLIST *mlist)
Definition ecat63ml.c:31
int mat_numcod(int frame, int plane, int gate, int data, int bed)
Definition ecat63ml.c:242
void mat_numdoc(int matnum, Matval *matval)
Definition ecat63ml.c:254
void ecat63PrintMainheader(ECAT63_mainheader *h, FILE *fp)
Definition ecat63p.c:16
int ecat63ReadMainheader(FILE *fp, ECAT63_mainheader *h)
Definition ecat63r.c:25
int ecat63WriteScan(FILE *fp, int matnum, ECAT63_scanheader *h, void *data)
Definition ecat63w.c:461
int ecat63WriteImage(FILE *fp, int matnum, ECAT63_imageheader *h, void *data)
Definition ecat63w.c:410
FILE * ecat63Create(const char *fname, ECAT63_mainheader *h)
Definition ecat63w.c:365
int ecat63WriteMainheader(FILE *fp, ECAT63_mainheader *h)
Definition ecat63w.c:24
int IMG_TEST
Definition img.c:6
void imgInfo(IMG *image)
Definition img.c:359
char * imgIsotope(IMG *img)
Definition imgdecayc.c:76
int imgUnitToEcat6(IMG *img)
Definition imgunits.c:175
#define IMG_TYPE_RAW
#define IMG_STATUS_OCCUPIED
#define ECAT63_SYSTEM_TYPE_DEFAULT
#define IMG_DC_CORRECTED
#define RAW_DATA
#define IMAGE_DATA
#define VAX_I2
#define IMG_TYPE_IMAGE
size_t strlcpy(char *dst, const char *src, size_t dstsize)
Definition strext.c:245
int temp_roundf(float e)
Definition petc99.c:20
short int num_dimensions
char radiopharmaceutical[32]
short int scan_start_month
short int scan_start_second
short int num_bed_pos
short int calibration_units
short int scan_start_year
char study_description[32]
short int scan_start_day
short int scan_start_minute
char user_process_code[10]
short int scan_start_hour
short int system_type
short int dimension_2
short int dimension_1
short int frame_duration_sec
float sizex
unsigned short int dimx
char type
char patientName[32]
char studyDescription[32]
float sampleDistance
float **** m
char decayCorrection
float transaxialFOV
char status
time_t scanStart
float * prompts
char userProcessCode[11]
unsigned short int dimt
int _dataType
int * planeNumber
char patientID[16]
int scanner
float * start
unsigned short int dimz
unsigned short int dimy
float * end
float zoom
char radiopharmaceutical[32]
float * decayCorrFactor
float isotopeHalflife
char studyNr[MAX_STUDYNR_LEN+1]
float * randoms
float axialFOV
float sizez
MatDir * matdir
int matnum
int frame
int plane

◆ ecat63ReadAllToImg()

int ecat63ReadAllToImg ( const char * fname,
IMG * img )

Read all matrices in ECAT file to memory.

Sinograms are dead-time corrected.

See also
imgInit
Returns
0 if ok, 1 invalid input, 3 failed to open file for reading, 4 failed to read main header, 5 failed to read matrix list, 6 matrix not found, 7 variable matrix sizes, 8 failed to read matrix sub header, 9 failed to allocate memory for data, 10 failed to read sub header, 11 failed to read matrix data.
Parameters
fnameName of the input ECAT 6.3 file.
imgData structure in which the file is read; must be initialized before using this function.

Definition at line 22 of file img_e63.c.

28 {
29 int m, ret, blkNr=0, dim_x=0, dim_y=0;
30 int frameNr, planeNr, del_nr=0;
31 int frame, plane, prev_frame, prev_plane, seqplane;
32 FILE *fp;
33 ECAT63_mainheader main_header;
34 MATRIXLIST mlist;
35 MatDir tmpMatdir;
36 Matval matval;
37 ECAT63_imageheader image_header;
38 ECAT63_scanheader scan_header;
39 ECAT63_attnheader attn_header;
40 ECAT63_normheader norm_header;
41 float scale=1.0;
42 short int *sptr;
43 char *mdata=NULL, *mptr;
44 /*int *iptr;*/
45
46
47 if(IMG_TEST) printf("ecat63ReadAllToImg(%s, *img)\n", fname);
48 /* Check the arguments */
49 if(fname==NULL || img==NULL || img->status!=IMG_STATUS_INITIALIZED) {
50 if(img) imgSetStatus(img, STATUS_FAULT);
51 return 1;
52 }
53
54 /* Open input CTI file for read */
55 if((fp=fopen(fname, "rb")) == NULL) {
56 imgSetStatus(img, STATUS_NOFILE);
57 return 3;
58 }
59
60 /* Read main header */
61 if((ret=ecat63ReadMainheader(fp, &main_header))) {
62 imgSetStatus(img, STATUS_NOMAINHEADER);
63 return 4;
64 }
65 if(IMG_TEST>5) ecat63PrintMainheader(&main_header, stdout);
66
67 /* Read matrix list and nr */
68 ecat63InitMatlist(&mlist);
69 ret=ecat63ReadMatlist(fp, &mlist, ECAT63_TEST-1);
70 if(ret) {
71 imgSetStatus(img, STATUS_NOMATLIST);
72 return 5;
73 }
74 if(mlist.matrixNr<=0) {
75 strcpy(ecat63errmsg, "matrix list is empty");
76 return 5;
77 }
78 /* Sort matrix list */
79 for(int i=0; i<mlist.matrixNr-1; i++)
80 for(int j=i+1; j<mlist.matrixNr; j++) {
81 if(mlist.matdir[i].matnum>mlist.matdir[j].matnum) {
82 tmpMatdir=mlist.matdir[i];
83 mlist.matdir[i]=mlist.matdir[j]; mlist.matdir[j]=tmpMatdir;
84 }
85 }
86 if(IMG_TEST>6) {
87 printf("matrix list after sorting:\n");
88 ecat63PrintMatlist(&mlist);
89 }
90
91 /* Trim extra frames */
92 if(main_header.num_frames>0) {
93 del_nr=ecat63DeleteLateFrames(&mlist, main_header.num_frames);
94 if(IMG_TEST && del_nr>0)
95 printf(" %d entries in matrix list will not be used.\n", del_nr);
96 }
97 /* Calculate the number of planes, frames and gates */
98 /* Check that all planes have equal nr of frames (gates) */
99 /* and check that frames (gates) are consequentially numbered */
100 prev_plane=plane=-1; prev_frame=frame=-1; frameNr=planeNr=0; ret=0;
101 for(int m=0; m<mlist.matrixNr; m++) if(mlist.matdir[m].matstat==1) {
102 /* get frame and plane */
103 mat_numdoc(mlist.matdir[m].matnum, &matval);
104 plane=matval.plane;
105 if(main_header.num_frames>=main_header.num_gates) frame=matval.frame;
106 else frame=matval.gate;
107 if(IMG_TEST>2) printf("frame=%d plane=%d\n", frame, plane);
108 /* did plane change? */
109 if(plane!=prev_plane) {
110 frameNr=1; planeNr++;
111 } else {
112 frameNr++;
113 /* In each plane, frame (gate) numbers must be continuous */
114 if(prev_frame>0 && frame!=prev_frame+1) {ret=1; break;}
115 }
116 prev_plane=plane; prev_frame=frame;
117 /* Calculate and check the size of one data matrix */
118 if(m==0) {
119 blkNr=mlist.matdir[m].endblk-mlist.matdir[m].strtblk;
120 } else if(blkNr!=mlist.matdir[m].endblk-mlist.matdir[m].strtblk) {
121 ret=2; break;
122 }
123 } /* next matrix */
124 if(IMG_TEST) printf("frameNr=%d planeNr=%d\n", frameNr, planeNr);
125 if(ret==1 || (frameNr*planeNr != mlist.matrixNr-del_nr)) {
126 strcpy(ecat63errmsg, "missing matrix");
127 ecat63EmptyMatlist(&mlist); fclose(fp);
128 return(6); /* this number is used in imgRead() */
129 }
130 if(ret==2) {
131 strcpy(ecat63errmsg, "matrix sizes are different");
132 ecat63EmptyMatlist(&mlist); fclose(fp); return(7);
133 }
134 /* Read x,y-dimensions from 1st matrix subheader */
135 m=0; if(main_header.file_type==IMAGE_DATA) {
136 ret=ecat63ReadImageheader(fp, mlist.matdir[m].strtblk, &image_header, IMG_TEST-4, NULL);
137 dim_x=image_header.dimension_1; dim_y=image_header.dimension_2;
138 } else if(main_header.file_type==RAW_DATA) {
139 ret=ecat63ReadScanheader(fp, mlist.matdir[m].strtblk, &scan_header, IMG_TEST-4, NULL);
140 dim_x=scan_header.dimension_1; dim_y=scan_header.dimension_2;
141 } else if(main_header.file_type==ATTN_DATA) {
142 ret=ecat63ReadAttnheader(fp, mlist.matdir[m].strtblk, &attn_header, IMG_TEST-10, NULL);
143 dim_x=attn_header.dimension_1; dim_y=attn_header.dimension_2;
144 } else if(main_header.file_type==NORM_DATA) {
145 ret=ecat63ReadNormheader(fp, mlist.matdir[m].strtblk, &norm_header, IMG_TEST-10, NULL);
146 dim_x=norm_header.dimension_1; dim_y=norm_header.dimension_2;
147 }
148 if(ret) {
149 sprintf(ecat63errmsg, "cannot read matrix %u subheader",
150 mlist.matdir[m].matnum);
151 ecat63EmptyMatlist(&mlist); fclose(fp); return(8);
152 }
153
154 /* Allocate memory for ECAT data matrix */
155 if(IMG_TEST>1) printf("allocating memory for %d blocks\n", blkNr);
156 mdata=(char*)malloc((size_t)blkNr*MatBLKSIZE); if(mdata==NULL) {
157 strcpy(ecat63errmsg, "out of memory");
158 fclose(fp); ecat63EmptyMatlist(&mlist); return(8);
159 }
160 /* Allocate memory for whole img data */
161 ret=imgAllocate(img, planeNr, dim_y, dim_x, frameNr);
162 if(ret) {
163 sprintf(ecat63errmsg, "out of memory (%d)", ret);
164 fclose(fp); ecat63EmptyMatlist(&mlist); free(mdata); return(9);
165 }
166
167 /* Fill img info with ECAT main and sub header information */
168 img->scanner=main_header.system_type;
169 img->unit=main_header.calibration_units; /* this continues below */
170 strlcpy(img->radiopharmaceutical, main_header.radiopharmaceutical, 32);
171 img->isotopeHalflife=main_header.isotope_halflife;
172 img->scanStart=ecat63Scanstarttime(&main_header);
173 if(img->scanStart==-1) {
174 img->scanStart=0;
175 if(IMG_TEST>0) printf("invalid scan_start_time in mainheader\n");
176 }
177 img->axialFOV=10.*main_header.axial_fov;
178 img->transaxialFOV=10.*main_header.transaxial_fov;
179 strlcpy(img->studyNr, main_header.study_name, MAX_STUDYNR_LEN+1);
180 strlcpy(img->patientName, main_header.patient_name, 32);
181 strlcpy(img->patientID, main_header.patient_id, 16);
182 img->sizez=10.*main_header.plane_separation;
183 if(main_header.file_type==IMAGE_DATA) {
184 img->type=IMG_TYPE_IMAGE;
185 if(img->unit<1) img->unit=image_header.quant_units;
186 img->zoom=image_header.recon_scale;
187 if(image_header.decay_corr_fctr>1.0) img->decayCorrection=IMG_DC_CORRECTED;
189 img->sizex=img->sizey=10.*image_header.pixel_size;
190 if(img->sizez<10.*image_header.slice_width)
191 img->sizez=10.*image_header.slice_width;
192 img->xform[0]=NIFTI_XFORM_UNKNOWN; // qform
193 img->xform[1]=NIFTI_XFORM_SCANNER_ANAT; // sform
194 img->quatern[6]=img->sizex; img->quatern[9]=img->sizex;
195 img->quatern[11]=img->sizey; img->quatern[13]=img->sizey;
196 img->quatern[16]=img->sizez; img->quatern[17]=img->sizez;
197 } else if(main_header.file_type==RAW_DATA) {
198 img->type=IMG_TYPE_RAW;
200 } else if(main_header.file_type==ATTN_DATA) {
201 img->type=IMG_TYPE_ATTN;
203 } else if(main_header.file_type==NORM_DATA) {
204 img->type=IMG_TYPE_RAW;
206 }
207 strlcpy(img->studyDescription, main_header.study_description, 32);
208 strncpy(img->userProcessCode, main_header.user_process_code, 10);
209 img->userProcessCode[10]=(char)0;
210 /* If valid study number is found in user_process_code, then take it */
213
214 /* Set _fileFormat */
215 img->_fileFormat=IMG_E63;
216
217 /* Read one ECAT matrix at a time and put them to img */
218 prev_plane=plane=-1; seqplane=-1;
219 for(int m=0; m<mlist.matrixNr; m++) if(mlist.matdir[m].matstat==1) {
220 /* get plane */
221 mat_numdoc(mlist.matdir[m].matnum, &matval);
222 plane=matval.plane;
223 /* did plane change? */
224 if(plane!=prev_plane) {seqplane++; frame=1;} else frame++;
225 prev_plane=plane;
226 img->planeNumber[seqplane]=plane;
227 /* Read subheader */
228 if(main_header.file_type==IMAGE_DATA) {
229 ret=ecat63ReadImageheader(fp, mlist.matdir[m].strtblk, &image_header, IMG_TEST-10, NULL);
230 if(ret==0 && (dim_x!=image_header.dimension_1 || dim_y!=image_header.dimension_2)) ret=1;
231 scale=image_header.quant_scale;
232 if(image_header.ecat_calibration_fctr>0.0 &&
233 fabs(main_header.calibration_factor/image_header.ecat_calibration_fctr - 1.0)>0.0001)
234 scale*=image_header.ecat_calibration_fctr;
235 if(img->unit==0 && image_header.quant_units>0)
236 img->unit=image_header.quant_units;
237 img->_dataType=image_header.data_type;
238 img->start[frame-1]=image_header.frame_start_time/1000.;
239 img->end[frame-1]=img->start[frame-1]+image_header.frame_duration/1000.;
240 img->mid[frame-1]=0.5*(img->start[frame-1]+img->end[frame-1]);
241 if(image_header.decay_corr_fctr>1.0)
242 img->decayCorrFactor[frame-1]=image_header.decay_corr_fctr;
243 else
244 img->decayCorrFactor[frame-1]=0.0;
245 /* Dead-time correction is assumed to be included in scale factor */
246 } else if(main_header.file_type==RAW_DATA) {
247 ret=ecat63ReadScanheader(fp, mlist.matdir[m].strtblk, &scan_header, IMG_TEST-10, NULL);
248 if(ret==0 && (dim_x!=scan_header.dimension_1 || dim_y!=scan_header.dimension_2)) ret=1;
249 scale=scan_header.scale_factor;
250 /* Dead-time correction is done here */
251 if(scan_header.loss_correction_fctr>0.0)
252 scale*=scan_header.loss_correction_fctr;
253 img->_dataType=scan_header.data_type;
254 img->start[frame-1]=scan_header.frame_start_time/1000.;
255 img->end[frame-1]=img->start[frame-1]+scan_header.frame_duration/1000.;
256 img->mid[frame-1]=0.5*(img->start[frame-1]+img->end[frame-1]);
257 img->sampleDistance=10.0*scan_header.sample_distance;
258 img->prompts[frame-1]=(float)scan_header.prompts;
259 img->randoms[frame-1]=(float)scan_header.delayed;
260 } else if(main_header.file_type==ATTN_DATA) {
261 ret=ecat63ReadAttnheader(fp, mlist.matdir[m].strtblk, &attn_header, IMG_TEST-10, NULL);
262 if(ret==0 && (dim_x!=attn_header.dimension_1 || dim_y!=attn_header.dimension_2)) ret=1;
263 img->_dataType=attn_header.data_type;
264 scale=attn_header.scale_factor;
265 img->sampleDistance=10.0*attn_header.sample_distance;
266 } else if(main_header.file_type==NORM_DATA) {
267 ret=ecat63ReadNormheader(fp, mlist.matdir[m].strtblk, &norm_header, IMG_TEST-10, NULL);
268 if(ret==0 && (dim_x!=norm_header.dimension_1 || dim_y!=norm_header.dimension_2)) ret=1;
269 scale=norm_header.scale_factor;
270 img->_dataType=norm_header.data_type;
271 } else img->_dataType=-1;
272 if(ret) {
273 sprintf(ecat63errmsg, "cannot read matrix %u subheader", mlist.matdir[m].matnum);
274 ecat63EmptyMatlist(&mlist); free(mdata); fclose(fp); return(10);
275 }
276 if(main_header.calibration_factor>0.0) scale*=main_header.calibration_factor;
277 if(IMG_TEST>2) printf("scale=%e datatype=%d\n", scale, img->_dataType);
278 /* Read the pixel data */
279 ret=ecat63ReadMatdata(fp, mlist.matdir[m].strtblk+1,
280 mlist.matdir[m].endblk-mlist.matdir[m].strtblk,
281 mdata, img->_dataType);
282 if(ret) {
283 strcpy(ecat63errmsg, "cannot read matrix data");
284 ecat63EmptyMatlist(&mlist); free(mdata); fclose(fp); return(11);
285 }
286 mptr=mdata;
287 if(img->_dataType==BYTE_TYPE) {
288 for(int i=0; i<dim_y; i++) for(int j=0; j<dim_x; j++, mptr++) {
289 img->m[seqplane][i][j][frame-1]=scale*(float)(*mptr);
290 }
291 } else if(img->_dataType==VAX_I2 || img->_dataType==SUN_I2) {
292 for(int i=0; i<dim_y; i++) for(int j=0; j<dim_x; j++, mptr+=2) {
293 sptr=(short int*)mptr;
294 img->m[seqplane][i][j][frame-1]=scale*(float)(*sptr);
295 }
296 } else if(img->_dataType==VAX_I4 || img->_dataType==SUN_I4) {
297 for(int i=0; i<dim_y; i++) for(int j=0; j<dim_x; j++, mptr+=4) {
298 /*iptr=(int*)mptr;*/
299 img->m[seqplane][i][j][frame-1]=1.0; /*scale*(float)(*iptr);*/
300 }
301 } else if(img->_dataType==VAX_R4 || img->_dataType==IEEE_R4) {
302 for(int i=0; i<dim_y; i++) for(int j=0; j<dim_x; j++, mptr+=4) {
303 memcpy(&img->m[seqplane][i][j][frame-1], mptr, 4);
304 img->m[seqplane][i][j][frame-1]*=scale;
305 }
306 }
307 } /* next matrix */
308
309 /* Set the unit */
310 imgUnitFromEcat(img, img->unit);
311
312 /* Free memory and close file */
313 free(mdata);
314 ecat63EmptyMatlist(&mlist);
315 fclose(fp);
316
317 /* For saving, only 2-byte VAX data types are allowed */
318 if(img->_dataType==VAX_I4 || img->_dataType==VAX_R4) img->_dataType=VAX_I2;
319
320 return(0);
321}
int ECAT63_TEST
Definition ecat63h.c:6
int ecat63DeleteLateFrames(MATRIXLIST *ml, int frame_nr)
Definition ecat63ml.c:342
void ecat63PrintMatlist(MATRIXLIST *ml)
Definition ecat63ml.c:130
int ecat63ReadNormheader(FILE *fp, int blk, ECAT63_normheader *h, int verbose, char *errmsg)
Definition ecat63r.c:483
int ecat63ReadScanheader(FILE *fp, int blk, ECAT63_scanheader *h, int verbose, char *errmsg)
Definition ecat63r.c:377
int ecat63ReadMatdata(FILE *fp, int strtblk, int blkNr, char *data, int dtype)
Definition ecat63r.c:568
int ecat63ReadAttnheader(FILE *fp, int blk, ECAT63_attnheader *h, int verbose, char *errmsg)
Definition ecat63r.c:296
int ecat63ReadImageheader(FILE *fp, int blk, ECAT63_imageheader *h, int verbose, char *errmsg)
Definition ecat63r.c:187
time_t ecat63Scanstarttime(const ECAT63_mainheader *h)
Get calendar time from ECAT 6.3 main header.
Definition ecat63w.c:925
int imgAllocate(IMG *image, int planes, int rows, int columns, int frames)
Definition img.c:194
void imgSetStatus(IMG *img, int status_index)
Definition img.c:345
void imgUnitFromEcat(IMG *img, int ecat_unit)
Definition imgunits.c:95
#define IMG_TYPE_ATTN
#define NIFTI_XFORM_UNKNOWN
#define SUN_I4
#define ATTN_DATA
#define VAX_R4
#define IMG_DC_UNKNOWN
#define IEEE_R4
#define IMG_STATUS_INITIALIZED
#define IMG_E63
#define IMG_DC_NONCORRECTED
#define BYTE_TYPE
#define VAX_I4
#define NORM_DATA
#define NIFTI_XFORM_SCANNER_ANAT
#define MatBLKSIZE
#define SUN_I2
#define MAX_STUDYNR_LEN
Definition libtpcmisc.h:163
int studynr_validity_check(char *studynr)
Definition studynr.c:196
short int dimension_1
short int dimension_2
short int dimension_1
short int dimension_2
short int xform[2]
char unit
int _fileFormat
float sizey
float quatern[18]
float * mid
int matstat
int endblk
int strtblk

Referenced by imgRead().

◆ ecat63ReadPlaneToImg()

int ecat63ReadPlaneToImg ( const char * fname,
IMG * img )

Reads one CTI ECAT 6.3 plane (all frames or gates) at a time to memory. Img data must be initialized before this procedure. Existing img->_dataType is not changed. If img data structure is empty, reads the first plane. If img data structure contains data, reads the next plane. Any existing data in img is cleared and replaced by the new plane.

Parameters
fnamename of the input ECAT 6.3 file
imgdata structure in which the file is read
Returns
0 if ok, 1 next plane was requested but not found any more, 2 invalid input data, 3 failed to open file, 4 failed to read main header, 5 failed to read matrix list, 6 invalid matrix data, 7 failed to read matrix sub header, 8 failed to allocate memory, 9 failed to read matrix data

Definition at line 550 of file img_e63.c.

550 {
551 int m, ret, blkNr=0, dim_x=0, dim_y=0, del_nr=0;
552 int frameNr, planeNr, datatype=0, firstm=0, isFirst=0;
553 int frame, plane, prev_frame, prev_plane, prev_frameNr=0;
554 FILE *fp;
555 ECAT63_mainheader main_header;
556 MATRIXLIST mlist;
557 MatDir tmpMatdir;
558 Matval matval;
559 ECAT63_imageheader image_header;
560 ECAT63_scanheader scan_header;
561 ECAT63_attnheader attn_header;
562 ECAT63_normheader norm_header;
563 float scale=1.0;
564 short int *sptr;
565 char *mdata=NULL, *mptr;
566 int *iptr;
567
568
569 if(IMG_TEST) printf("ecat63ReadPlaneToImg(%s, *img)\n", fname);
570 /* Check the arguments */
571 if(fname==NULL) {strcpy(ecat63errmsg, "invalid ECAT filename"); return(2);}
572 if(img==NULL || img->status==IMG_STATUS_UNINITIALIZED) {
573 strcpy(ecat63errmsg, "image data not initialized"); return(2);}
574
575 /* Open input CTI file for read */
576 if((fp=fopen(fname, "rb")) == NULL) {
577 strcpy(ecat63errmsg, "cannot open ECAT file"); return(3);}
578
579 /* Read main header */
580 if((ret=ecat63ReadMainheader(fp, &main_header))) {
581 sprintf(ecat63errmsg, "cannot read main header (%d)", ret);
582 fclose(fp); return(4);
583 }
584 if(IMG_TEST) ecat63PrintMainheader(&main_header, stdout);
585
586 /* Read matrix list and nr */
587 ecat63InitMatlist(&mlist);
588 ret=ecat63ReadMatlist(fp, &mlist, IMG_TEST);
589 if(ret) {
590 sprintf(ecat63errmsg, "cannot read matrix list (%d)", ret);
591 fclose(fp); return(5);
592 }
593 if(mlist.matrixNr<=0) {
594 strcpy(ecat63errmsg, "matrix list is empty"); fclose(fp); return(5);}
595 /* Sort matrix list */
596 for(int i=0; i<mlist.matrixNr-1; i++) for(int j=i+1; j<mlist.matrixNr; j++) {
597 if(mlist.matdir[i].matnum>mlist.matdir[j].matnum) {
598 tmpMatdir=mlist.matdir[i];
599 mlist.matdir[i]=mlist.matdir[j]; mlist.matdir[j]=tmpMatdir;
600 }
601 }
602 /* Trim extra frames */
603 if(main_header.num_frames>0) {
604 del_nr=ecat63DeleteLateFrames(&mlist, main_header.num_frames);
605 if(IMG_TEST && del_nr>0)
606 printf(" %d entries in matrix list will not be used.\n", del_nr);
607 }
608
609 /* Decide the plane to read */
610 if(img->status==IMG_STATUS_OCCUPIED) {
611 /* img contains data */
612 isFirst=0; prev_frameNr=img->dimt;
613 /* get the current plane in there */
614 prev_plane=img->planeNumber[img->dimz-1];
615 /* Clear all data in img: not here but only after finding new data */
616 } else {
617 /* img does not contain any data */
618 isFirst=1;
619 /* set "current plane" to -1 */
620 prev_plane=-1;
621 }
622 /* from sorted matrix list, find first plane larger than the current one */
623 plane=-1;
624 for(int m=0; m<mlist.matrixNr; m++) if(mlist.matdir[m].matstat==1) {
625 mat_numdoc(mlist.matdir[m].matnum, &matval);
626 if(matval.plane>prev_plane) {plane=matval.plane; break;}
627 }
628 /* If not found, return an error or 1 if this was not the first one */
629 if(plane<0) {
630 fclose(fp); ecat63EmptyMatlist(&mlist);
631 if(isFirst) {
632 sprintf(ecat63errmsg, "ECAT file contains no matrices");
633 return(7);
634 } else {
635 sprintf(ecat63errmsg, "ECAT file contains no more planes");
636 if(IMG_TEST) printf("%s\n", ecat63errmsg);
637 return(1);
638 }
639 }
640 if(IMG_TEST) printf("Next plane to read: %d\n", plane);
641 /* Clear all data in img */
642 imgEmpty(img);
643
644 /* Calculate the number of frames and gates */
645 prev_frame=frame=-1; frameNr=0; ret=0; planeNr=1;
646 for(int m=0; m<mlist.matrixNr; m++) if(mlist.matdir[m].matstat==1) {
647 /* get frame and plane; work only with current plane */
648 mat_numdoc(mlist.matdir[m].matnum, &matval);
649 if(matval.plane<plane) continue; else if(matval.plane>plane) break;
650 if(main_header.num_frames>=main_header.num_gates) frame=matval.frame;
651 else frame=matval.gate;
652 frameNr++;
653 /* frame (gate) numbers must be continuous */
654 if(prev_frame>0 && frame!=prev_frame+1) {ret=1; break;}
655 prev_frame=frame;
656 /* Calculate and check the size of one data matrix */
657 if(frameNr==1) {
658 firstm=m;
659 blkNr=mlist.matdir[m].endblk-mlist.matdir[m].strtblk;
660 } else if(blkNr!=mlist.matdir[m].endblk-mlist.matdir[m].strtblk) {
661 ret=2; break;
662 }
663 prev_frame=frame;
664 } /* next matrix */
665 if(ret==1) {
666 strcpy(ecat63errmsg, "missing matrix");
667 ecat63EmptyMatlist(&mlist); fclose(fp); return(6);
668 }
669 if(ret==2) {
670 strcpy(ecat63errmsg, "matrix sizes are different");
671 ecat63EmptyMatlist(&mlist); fclose(fp); return(6);
672 }
673 /* Check that frameNr is equal to the dimt in IMG */
674 if(!isFirst && frameNr!=prev_frameNr) {
675 strcpy(ecat63errmsg, "different frame nr in different planes");
676 ecat63EmptyMatlist(&mlist); fclose(fp); return(6);
677 }
678
679 /* Read x,y-dimensions from 1st matrix subheader on current plane */
680 m=firstm; if(main_header.file_type==IMAGE_DATA) {
681 ret=ecat63ReadImageheader(fp, mlist.matdir[m].strtblk, &image_header, IMG_TEST-10, NULL);
682 dim_x=image_header.dimension_1; dim_y=image_header.dimension_2;
683 } else if(main_header.file_type==RAW_DATA) {
684 ret=ecat63ReadScanheader(fp, mlist.matdir[m].strtblk, &scan_header, IMG_TEST-10, NULL);
685 dim_x=scan_header.dimension_1; dim_y=scan_header.dimension_2;
686 } else if(main_header.file_type==ATTN_DATA) {
687 ret=ecat63ReadAttnheader(fp, mlist.matdir[m].strtblk, &attn_header, IMG_TEST-10, NULL);
688 dim_x=attn_header.dimension_1; dim_y=attn_header.dimension_2;
689 } else if(main_header.file_type==NORM_DATA) {
690 ret=ecat63ReadNormheader(fp, mlist.matdir[m].strtblk, &norm_header, IMG_TEST-10, NULL);
691 dim_x=norm_header.dimension_1; dim_y=norm_header.dimension_2;
692 }
693 if(ret) {
694 sprintf(ecat63errmsg, "cannot read matrix %u subheader",
695 mlist.matdir[m].matnum);
696 ecat63EmptyMatlist(&mlist); fclose(fp); return(7);
697 }
698
699 /* Allocate memory for ECAT data matrix */
700 if(IMG_TEST) printf("allocating memory for %d blocks\n", blkNr);
701 mdata=(char*)malloc((size_t)blkNr*MatBLKSIZE); if(mdata==NULL) {
702 strcpy(ecat63errmsg, "out of memory");
703 fclose(fp); ecat63EmptyMatlist(&mlist); return(8);
704 }
705 /* Allocate memory for whole img data */
706 ret=imgAllocate(img, planeNr, dim_y, dim_x, frameNr);
707 if(ret) {
708 sprintf(ecat63errmsg, "out of memory (%d)", ret);
709 fclose(fp); ecat63EmptyMatlist(&mlist); free(mdata); return(8);
710 }
711
712 /* Fill img info with ECAT main and sub header information */
713 img->scanner=main_header.system_type;
714 img->unit=main_header.calibration_units; /* this continues below */
715 strlcpy(img->radiopharmaceutical, main_header.radiopharmaceutical, 32);
716 img->isotopeHalflife=main_header.isotope_halflife;
717 img->scanStart=ecat63Scanstarttime(&main_header);
718 if(img->scanStart==-1) {
719 img->scanStart=0;
720 if(IMG_TEST>0) printf("invalid scan_start_time in mainheader\n");
721 }
722 img->axialFOV=10.*main_header.axial_fov;
723 img->transaxialFOV=10.*main_header.transaxial_fov;
724 strlcpy(img->studyNr, main_header.study_name, MAX_STUDYNR_LEN);
725 strlcpy(img->patientName, main_header.patient_name, 32);
726 strlcpy(img->patientID, main_header.patient_id, 16);
727 img->sizez=10.*main_header.plane_separation;
728 if(main_header.file_type==IMAGE_DATA) {
729 img->type=IMG_TYPE_IMAGE;
730 if(img->unit<1) img->unit=image_header.quant_units;
731 img->zoom=image_header.recon_scale;
732 if(image_header.decay_corr_fctr>1.0) img->decayCorrection=IMG_DC_CORRECTED;
734 img->sizex=img->sizey=10.*image_header.pixel_size;
735 if(img->sizez<10.*image_header.slice_width)
736 img->sizez=10.*image_header.slice_width;
737 img->xform[0]=NIFTI_XFORM_UNKNOWN; // qform
738 img->xform[1]=NIFTI_XFORM_SCANNER_ANAT; // sform
739 img->quatern[6]=img->sizex; img->quatern[9]=img->sizex;
740 img->quatern[11]=img->sizey; img->quatern[13]=img->sizey;
741 img->quatern[16]=img->sizez; img->quatern[17]=img->sizez;
742 } else if(main_header.file_type==RAW_DATA) {
743 img->type=IMG_TYPE_RAW;
745 } else if(main_header.file_type==ATTN_DATA) {
746 img->type=IMG_TYPE_ATTN;
748 } else if(main_header.file_type==NORM_DATA) {
749 img->type=IMG_TYPE_RAW;
751 }
752 img->planeNumber[0]=plane;
753 strlcpy(img->studyDescription, main_header.study_description, 32);
754 strlcpy(img->userProcessCode, main_header.user_process_code, 10);
755 //img->userProcessCode[10]=(char)0;
756 /* If valid study number is found in user_process_code, then take it */
759 /* Set _fileFormat */
760 img->_fileFormat=IMG_E63;
761
762 /* Read one ECAT matrix at a time and put them to img */
763 for(m=firstm, frame=1; m<mlist.matrixNr; m++) if(mlist.matdir[m].matstat==1) {
764 /* get plane */
765 mat_numdoc(mlist.matdir[m].matnum, &matval);
766 if(matval.plane>plane) break; /* Quit when current plane is read */
767 /* Read subheader */
768 if(main_header.file_type==IMAGE_DATA) {
769 ret=ecat63ReadImageheader(fp, mlist.matdir[m].strtblk, &image_header, IMG_TEST-10, NULL);
770 if(ret==0 && (dim_x!=image_header.dimension_1 || dim_y!=image_header.dimension_2)) ret=1;
771 scale=image_header.quant_scale;
772 if(image_header.ecat_calibration_fctr>0.0 &&
773 fabs(main_header.calibration_factor/image_header.ecat_calibration_fctr - 1.0)>0.0001)
774 scale*=image_header.ecat_calibration_fctr;
775 if(img->unit==0 && image_header.quant_units>0)
776 img->unit=image_header.quant_units;
777 datatype=image_header.data_type;
778 img->start[frame-1]=image_header.frame_start_time/1000.;
779 img->end[frame-1]=img->start[frame-1]+image_header.frame_duration/1000.;
780 img->mid[frame-1]=0.5*(img->start[frame-1]+img->end[frame-1]);
781 if(image_header.decay_corr_fctr>1.0)
782 img->decayCorrFactor[frame-1]=image_header.decay_corr_fctr;
783 else
784 img->decayCorrFactor[frame-1]=0.0;
785 } else if(main_header.file_type==RAW_DATA) {
786 ret=ecat63ReadScanheader(fp, mlist.matdir[m].strtblk, &scan_header, IMG_TEST-10, NULL);
787 if(ret==0 && (dim_x!=scan_header.dimension_1 ||dim_y!=scan_header.dimension_2)) ret=1;
788 scale=scan_header.scale_factor;
789 if(scan_header.loss_correction_fctr>0.0) // dead-time correction
790 scale*=scan_header.loss_correction_fctr;
791 datatype=scan_header.data_type;
792 img->start[frame-1]=scan_header.frame_start_time/1000.;
793 img->end[frame-1]=img->start[frame-1]+scan_header.frame_duration/1000.;
794 img->mid[frame-1]=0.5*(img->start[frame-1]+img->end[frame-1]);
795 img->sampleDistance=10.0*scan_header.sample_distance;
796 img->prompts[frame-1]=(float)scan_header.prompts;
797 img->randoms[frame-1]=(float)scan_header.delayed;
798 } else if(main_header.file_type==ATTN_DATA) {
799 ret=ecat63ReadAttnheader(fp, mlist.matdir[m].strtblk, &attn_header, IMG_TEST-10, NULL);
800 if(ret==0 && (dim_x!=attn_header.dimension_1 ||dim_y!=attn_header.dimension_2)) ret=1;
801 img->sampleDistance=10.0*attn_header.sample_distance;
802 scale=attn_header.scale_factor;
803 datatype=attn_header.data_type;
804 } else if(main_header.file_type==NORM_DATA) {
805 ret=ecat63ReadNormheader(fp, mlist.matdir[m].strtblk, &norm_header, IMG_TEST-10, NULL);
806 if(ret==0 && (dim_x!=norm_header.dimension_1 ||dim_y!=norm_header.dimension_2)) ret=1;
807 scale=norm_header.scale_factor;
808 datatype=norm_header.data_type;
809 } else datatype=-1;
810 if(ret) {
811 sprintf(ecat63errmsg, "cannot read matrix %u subheader",
812 mlist.matdir[m].matnum);
813 ecat63EmptyMatlist(&mlist); free(mdata); fclose(fp); return(7);
814 }
815 if(main_header.calibration_factor>0.0)
816 scale*=main_header.calibration_factor;
817 if(IMG_TEST>2) printf("scale=%e datatype=%d\n", scale, datatype);
818 /* Read the pixel data */
819 ret=ecat63ReadMatdata(fp, mlist.matdir[m].strtblk+1,
820 mlist.matdir[m].endblk-mlist.matdir[m].strtblk,
821 mdata, datatype);
822 if(ret) {
823 strcpy(ecat63errmsg, "cannot read matrix data");
824 ecat63EmptyMatlist(&mlist); free(mdata); fclose(fp); return(9);
825 }
826 mptr=mdata;
827 if(datatype==BYTE_TYPE) {
828 for(int i=0; i<dim_y; i++) for(int j=0; j<dim_x; j++, mptr++) {
829 img->m[0][i][j][frame-1]=scale*(float)(*mptr);
830 }
831 } else if(datatype==VAX_I2 || datatype==SUN_I2) {
832 for(int i=0; i<dim_y; i++) for(int j=0; j<dim_x; j++, mptr+=2) {
833 sptr=(short int*)mptr;
834 img->m[0][i][j][frame-1]=scale*(float)(*sptr);
835 }
836 } else if(datatype==VAX_I4 || datatype==SUN_I4) {
837 for(int i=0; i<dim_y; i++) for(int j=0; j<dim_x; j++, mptr+=4) {
838 iptr=(int*)mptr;
839 img->m[0][i][j][frame-1]=scale*(float)(*iptr);
840 }
841 } else if(datatype==VAX_R4 || datatype==IEEE_R4) {
842 for(int i=0; i<dim_y; i++) for(int j=0; j<dim_x; j++, mptr+=4) {
843 memcpy(&img->m[0][i][j][frame-1], mptr, 4);
844 img->m[0][i][j][frame-1]*=scale;
845 }
846 }
847 frame++;
848 } /* next matrix (frame) */
849 /* Set the unit */
850 imgUnitFromEcat(img, img->unit);
851
852 /* Free memory and close file */
853 free(mdata);
854 ecat63EmptyMatlist(&mlist);
855 fclose(fp);
856
857 /* Set _dataType if it has not been set before */
858 if(img->_dataType<1) img->_dataType=datatype;
859 /* For saving, only 2-byte VAX data types are allowed */
860 if(img->_dataType==VAX_I4 || img->_dataType==VAX_R4) img->_dataType=VAX_I2;
861
862 return(0);
863}
void imgEmpty(IMG *image)
Definition img.c:121
#define IMG_STATUS_UNINITIALIZED

◆ ecat63WriteAllImg()

int ecat63WriteAllImg ( const char * fname,
IMG * img )

Write all matrices in memory to the ECAT file.

Parameters
fnamename of the output ECAT 6.3 file, If ECAT file exists, it is renamed as filename%
imgdata structure from which the data is written
Returns
0 if ok, 1 invalid data, 2 image status is not 'occupied', 3 failed to create file, 4 failed to allocate memory for data, 9 failed to write data

Definition at line 335 of file img_e63.c.

335 {
336 int frame, plane, m, ret=0;
337 float f, fmax, fmin, g, scale;
338 short int *sdata, *sptr, smin, smax;
339 FILE *fp;
340 ECAT63_mainheader main_header;
341 ECAT63_imageheader image_header;
342 ECAT63_scanheader scan_header;
343 struct tm tm;
344
345 if(IMG_TEST) printf("ecat63WriteAllImg(%s, *img)\n", fname);
346 /* Check the arguments */
347 if(fname==NULL) {strcpy(ecat63errmsg, "invalid ECAT filename"); return(1);}
348 if(img==NULL || img->status!=IMG_STATUS_OCCUPIED) {
349 strcpy(ecat63errmsg, "image data is empty"); return(2);}
350 if(img->_dataType<1) img->_dataType=VAX_I2;
351
352 /* Check image size */
353 if(img->dimt>511 || img->dimz>255 || img->dimx>1024 || img->dimy>1024) {
354 strcpy(ecat63errmsg, "too large matrix size"); return(1);}
355
356 /* Initiate headers */
357 memset(&main_header, 0, sizeof(ECAT63_mainheader));
358 memset(&image_header,0, sizeof(ECAT63_imageheader));
359 memset(&scan_header, 0, sizeof(ECAT63_scanheader));
360
361 /* Make a main header */
362 main_header.sw_version=2;
363 main_header.num_planes=img->dimz;
364 main_header.num_frames=img->dimt;
365 main_header.num_gates=1;
366 main_header.num_bed_pos=1;
367 if(img->type==IMG_TYPE_IMAGE) main_header.file_type=IMAGE_DATA;
368 else if(img->type==IMG_TYPE_RAW) main_header.file_type=RAW_DATA;
369 else {strcpy(ecat63errmsg, "invalid filetype"); return(1);}
370 main_header.data_type=img->_dataType;
371 if(img->scanner>0) main_header.system_type=img->scanner;
373 main_header.calibration_factor=1.0;
374 main_header.axial_fov=img->axialFOV/10.0;
375 main_header.transaxial_fov=img->transaxialFOV/10.0;
376 main_header.plane_separation=img->sizez/10.0;
377 main_header.calibration_units=imgUnitToEcat6(img);
378 strncpy(main_header.radiopharmaceutical, img->radiopharmaceutical, 32);
379 if(gmtime_r((time_t*)&img->scanStart, &tm)!=NULL) {
380 main_header.scan_start_year=tm.tm_year+1900;
381 main_header.scan_start_month=tm.tm_mon+1;
382 main_header.scan_start_day=tm.tm_mday;
383 main_header.scan_start_hour=tm.tm_hour;
384 main_header.scan_start_minute=tm.tm_min;
385 main_header.scan_start_second=tm.tm_sec;
386 if(IMG_TEST>2) {
387 printf(" img->scanStart := %ld\n", (long int)img->scanStart);
388 printf(" -> tm_year := %d\n", tm.tm_year);
389 printf(" -> tm_hour := %d\n", tm.tm_hour);
390 }
391 } else {
392 main_header.scan_start_year=1900;
393 main_header.scan_start_month=1;
394 main_header.scan_start_day=1;
395 main_header.scan_start_hour=0;
396 main_header.scan_start_minute=0;
397 main_header.scan_start_second=0;
398 if(IMG_TEST>0) printf("invalid scan_start_time in IMG\n");
399 }
400 main_header.isotope_halflife=img->isotopeHalflife;
401 strcpy(main_header.isotope_code, imgIsotope(img));
402 strlcpy(main_header.study_name, img->studyNr, 12);
403 strcpy(main_header.patient_name, img->patientName);
404 strcpy(main_header.patient_id, img->patientID);
405 strlcpy(main_header.user_process_code, img->userProcessCode, 10);
406 strncpy(main_header.study_description, img->studyDescription, 32);
407 if(IMG_TEST) ecat63PrintMainheader(&main_header, stdout);
408
409 /* Allocate memory for matrix data */
410 sdata=(short int*)malloc((size_t)img->dimx*img->dimy*sizeof(short int));
411 if(sdata==NULL) {strcpy(ecat63errmsg, "out of memory"); return(4);}
412
413 /* Open output ECAT file */
414 fp=ecat63Create(fname, &main_header);
415 if(fp==NULL) {strcpy(ecat63errmsg, "cannot write ECAT file"); return(3);}
416
417 /* Set the subheader contents, as far as possible */
418 switch(main_header.file_type) {
419 case RAW_DATA:
420 scan_header.data_type=main_header.data_type;
421 scan_header.dimension_1=img->dimx;
422 scan_header.dimension_2=img->dimy;
423 scan_header.frame_duration_sec=1;
424 scan_header.scale_factor=1.0;
425 scan_header.frame_start_time=0;
426 scan_header.frame_duration=1000;
427 /* Deadtime correction was done when reading, set to 1.0 to prevent
428 second correction when reading again. */
429 scan_header.loss_correction_fctr=1.0;
430 /*if(IMG_TEST) ecat63PrintScanheader(&scan_header);*/
431 break;
432 case IMAGE_DATA:
433 image_header.data_type=main_header.data_type;
434 image_header.num_dimensions=2;
435 image_header.dimension_1=img->dimx;
436 image_header.dimension_2=img->dimy;
437 image_header.recon_scale=img->zoom;
438 image_header.quant_scale=1.0;
439 image_header.slice_width=img->sizez/10.;
440 image_header.pixel_size=img->sizex/10.;
441 image_header.frame_start_time=0;
442 image_header.frame_duration=1000;
443 image_header.plane_eff_corr_fctr=1.0;
444 image_header.decay_corr_fctr=1.0;
445 image_header.loss_corr_fctr=1.0;
446 image_header.ecat_calibration_fctr=1.0;
447 image_header.well_counter_cal_fctr=1.0;
448 image_header.quant_units=main_header.calibration_units;
449 /*if(IMG_TEST) ecat63PrintImageheader(&image_header);*/
450 break;
451 }
452
453 /* Write one matrix at a time */
454 int n=0;
455 for(plane=1; plane<=img->dimz;plane++) for(frame=1;frame<=img->dimt;frame++) {
456 /* Scale data to short ints */
457 /* Search min and max */
458 fmin=fmax=f=img->m[plane-1][0][0][frame-1];
459 for(int i=0; i<img->dimy; i++) for(int j=0; j<img->dimx; j++) {
460 f=img->m[plane-1][i][j][frame-1];
461 if(f>fmax) fmax=f; else if(f<fmin) fmin=f;
462 }
463 /* Calculate scaling factor */
464 if(fabs(fmin)>fabs(fmax)) g=fabs(fmin); else g=fabs(fmax);
465 if(!isfinite(g)) {
466 sprintf(ecat63errmsg, "invalid pixel values on pl%02d fr%02d (%d).", plane, frame, ret);
467 fclose(fp); remove(fname); free(sdata);
468 return(1);
469 }
470 if(g>0.0) scale=32766./g; else scale=1.0;
471 if(!isfinite(scale)) scale=1.0;
472 // printf("fmin=%e fmax=%e g=%e scale=%e\n", fmin, fmax, g, scale);
473 /* Scale matrix data to shorts */
474 sptr=sdata;
475 for(int i=0; i<img->dimy; i++) for(int j=0; j<img->dimx; j++) {
476 *sptr=(short int)temp_roundf(scale*img->m[plane-1][i][j][frame-1]);
477 sptr++;
478 }
479 /* Calculate and set subheader min&max and scale */
480 smin=(short int)temp_roundf(scale*fmin);
481 smax=(short int)temp_roundf(scale*fmax);
482 if(main_header.file_type==RAW_DATA) {
483 scan_header.scan_min=smin; scan_header.scan_max=smax;
484 scan_header.scale_factor=1.0/scale;
485 scan_header.frame_start_time=(int)temp_roundf(1000.*img->start[frame-1]);
486 scan_header.frame_duration=
487 (int)temp_roundf(1000.*(img->end[frame-1]-img->start[frame-1]));
488 scan_header.sample_distance=(img->sampleDistance)/10.0;
489 scan_header.prompts=temp_roundf(img->prompts[frame-1]);
490 scan_header.delayed=temp_roundf(img->randoms[frame-1]);
491 } else if(main_header.file_type==IMAGE_DATA) {
492 image_header.image_min=smin; image_header.image_max=smax;
493 image_header.quant_scale=1.0/scale;
494 //if(!isfinite(image_header.quant_scale))
495 // printf("g=%g scale=%g quant_scale=%g\n", g, scale, image_header.quant_scale);
496 image_header.frame_start_time=(int)temp_roundf(1000.*img->start[frame-1]);
497 image_header.frame_duration=
498 (int)temp_roundf(1000.*(img->end[frame-1]-img->start[frame-1]));
499 /* Set decay correction factor */
501 image_header.decay_corr_fctr=img->decayCorrFactor[frame-1];
502 else
503 image_header.decay_corr_fctr=0.0;
504 }
505 /* Write matrix data */
506 m=mat_numcod(frame, img->planeNumber[plane-1], 1, 0, 0);
507 sptr=sdata;
508 if(IMG_TEST) printf(" writing matnum=%d\n", m);
509 if(main_header.file_type==RAW_DATA) {
510 if(IMG_TEST) ecat63PrintScanheader(&scan_header, stdout);
511 ret=ecat63WriteScan(fp, m, &scan_header, sptr);
512 } else if(main_header.file_type==IMAGE_DATA) {
513 if(IMG_TEST) ecat63PrintImageheader(&image_header, stdout);
514 ret=ecat63WriteImage(fp, m, &image_header, sptr);
515 }
516 if(ret) {
517 sprintf(ecat63errmsg, "cannot write data on pl%02d fr%02d (%d).",
518 plane, frame, ret);
519 fclose(fp); remove(fname); free(sdata);
520 return(9);
521 }
522 n++;
523 } /* next matrix */
524 if(IMG_TEST) printf(" %d matrices written in %s\n", n, fname);
525
526 /* Close file and free memory */
527 fclose(fp); free(sdata);
528
529 return(0);
530}
void ecat63PrintImageheader(ECAT63_imageheader *h, FILE *fp)
Definition ecat63p.c:94
void ecat63PrintScanheader(ECAT63_scanheader *h, FILE *fp)
Definition ecat63p.c:137

Referenced by imgWrite().

◆ imgEcat63Supported()

int imgEcat63Supported ( ECAT63_mainheader * h)

Check whether read functions in IMG library support this ECAT 6.3 file_type.

Parameters
hEcat 6.3 main header
Returns
1 if supported, 0 if not.

Definition at line 1155 of file img_e63.c.

1155 {
1156 if(h==NULL) return(0);
1157 if(h->file_type==IMAGE_DATA) return(1);
1158 if(h->file_type==RAW_DATA) return(1);
1159 if(h->file_type==ATTN_DATA) return(1);
1160 if(h->file_type==NORM_DATA) return(1);
1161 return(0);
1162}

Referenced by imgReadEcat63Header().

◆ imgGetEcat63Fileformat()

int imgGetEcat63Fileformat ( ECAT63_mainheader * h)

Return the IMG fileformat based on ECAT 6.3 file_type.

Parameters
hEcat 6.3 main header
Returns
IMG._fileFormat value.

Definition at line 1288 of file img_e63.c.

1288 {
1289 int fileFormat=IMG_UNKNOWN;
1290 switch(h->file_type) {
1291 case IMAGE_DATA:
1292 case RAW_DATA:
1293 case ATTN_DATA:
1294 case NORM_DATA:
1295 fileFormat=IMG_E63; break;
1296 default:
1297 fileFormat=IMG_UNKNOWN; break;
1298 }
1299 return fileFormat;
1300}
#define IMG_UNKNOWN

Referenced by imgReadEcat63Header().

◆ imgGetEcat63MHeader()

void imgGetEcat63MHeader ( IMG * img,
ECAT63_mainheader * h )

Copy ECAT 6.3 main header information into IMG

Parameters
imgtarget image structure
hsource Ecat 6.3 main header

Definition at line 1172 of file img_e63.c.

1173{
1174 if(IMG_TEST>0) printf("imgGetEcat63MHeader()\n");
1175 img->_dataType=h->data_type; /* again in subheaders*/
1176 /* For saving IMG data, only 2-byte VAX data types are allowed,
1177 so change it now */
1178 if(img->_dataType==VAX_I4 || img->_dataType==VAX_R4) img->_dataType=VAX_I2;
1179 img->scanner=h->system_type;
1183 if(img->scanStart==-1) {
1184 img->scanStart=0;
1185 if(IMG_TEST>0) printf("invalid scan_start_time in mainheader\n");
1186 }
1187 if(IMG_TEST>1) {
1188 char buf1[32], buf2[32];
1189 printf(" %s -> %s\n", ecat63ScanstarttimeInt(h, buf1),
1190 ctime_r_int(&img->scanStart, buf2));
1191 }
1192 img->axialFOV=10.0*h->axial_fov;
1193 img->transaxialFOV=10.0*h->transaxial_fov;
1195 strlcpy(img->patientName, h->patient_name, 32);
1196 strlcpy(img->patientID, h->patient_id, 16);
1197 img->sizez=10.0*h->plane_separation;
1198 switch(h->file_type) {
1199 case IMAGE_DATA: img->type=IMG_TYPE_IMAGE; break;
1200 case RAW_DATA:
1201 case ATTN_DATA:
1202 case NORM_DATA:
1203 default: img->type=IMG_TYPE_RAW;
1204 }
1206 strncpy(img->userProcessCode, h->user_process_code, 10);
1207 img->userProcessCode[10]=(char)0;
1208 /* If valid study number is found in user_process_code, then take it */
1209 if(!img->studyNr[0] && studynr_validity_check(img->userProcessCode))
1211 /* Set calibration unit (this may have to be read from subheader later) */
1213}
char * ctime_r_int(const time_t *t, char *buf)
Convert calendard time t into a null-terminated string of the form YYYY-MM-DD hh:mm:ss,...
Definition datetime.c:110
char * ecat63ScanstarttimeInt(const ECAT63_mainheader *h, char *buf)
Convert scan_start_time in ECAT 6.3 main header into a null-terminated string of the form YYYY-MM-DD ...
Definition ecat63p.c:391

Referenced by imgReadEcat63Header().

◆ imgReadEcat63FirstFrame()

int imgReadEcat63FirstFrame ( const char * fname,
IMG * img )

Read the first frame from an ECAT 6.3 file into IMG data structure.

Parameters
fnamename of file from which IMG contents will be read
imgpointer to the initiated but not preallocated IMG data
Returns
errstatus, which is STATUS_OK (0) when call was successful, and >0 in case of an error.

Definition at line 1457 of file img_e63.c.

1457 {
1458 int ret=0;
1459
1460 if(IMG_TEST) printf("\nimgReadEcat63FirstFrame(%s, *img)\n", fname);
1461 /* Check the input */
1462 if(img==NULL) return STATUS_FAULT;
1463 if(img->status!=IMG_STATUS_INITIALIZED) return STATUS_FAULT;
1464 imgSetStatus(img, STATUS_FAULT);
1465 if(fname==NULL) return STATUS_FAULT;
1466
1467 /* Read header information from file */
1468 ret=imgReadEcat63Header(fname, img); if(ret) return(ret);
1469 if(IMG_TEST>4) imgInfo(img);
1470
1471 /* Allocate memory for one frame */
1472 img->dimt=1;
1473 ret=imgAllocate(img, img->dimz, img->dimy, img->dimx, img->dimt);
1474 if(ret) return STATUS_NOMEMORY;
1475
1476 /* Read the first frame */
1477 ret=imgReadEcat63Frame(fname, 1, img, 0); if(ret) return(ret);
1478
1479 /* All went well */
1480 imgSetStatus(img, STATUS_OK);
1481 return STATUS_OK;
1482}
int imgReadEcat63Header(const char *fname, IMG *img)
Definition img_e63.c:1318
int imgReadEcat63Frame(const char *fname, int frame_to_read, IMG *img, int frame_index)
Definition img_e63.c:1503

◆ imgReadEcat63Frame()

int imgReadEcat63Frame ( const char * fname,
int frame_to_read,
IMG * img,
int frame_index )

Read a specified frame from an ECAT 6.3 file into preallocated IMG data structure.

IMG header is assumed to be filled correctly before calling this function, except for information concerning separate planes and this frame, which is filled here.

Parameters
fnamename of file from which IMG contents will be read
frame_to_readframe which will be read (1..frameNr)
imgpointer to the IMG data. Place for the frame must be preallocated
frame_indexIMG frame index (0..dimt-1) where data will be placed
Returns
errstatus, which is STATUS_OK (0) when call was successful, and >0 in case of an error. If frame does not exist, then and only then STATUS_NOMATRIX is returned.

Definition at line 1503 of file img_e63.c.

1505 {
1506 int ret=0, blkNr=0, frame, plane, seqplane=-1;
1507 int local_data_type=0;
1508 FILE *fp;
1509 ECAT63_mainheader main_header;
1510 MATRIXLIST mlist;
1511 Matval matval;
1512 ECAT63_imageheader image_header;
1513 ECAT63_scanheader scan_header;
1514 ECAT63_attnheader attn_header;
1515 ECAT63_normheader norm_header;
1516 float scale=1.0;
1517 short int *sptr;
1518 char *mdata=NULL, *mptr;
1519 /*int *iptr;*/
1520
1521
1522 if(IMG_TEST) printf("\nimgReadEcat63Frame(%s, %d, *img, %d)\n",
1523 fname, frame_to_read, frame_index);
1524
1525 /* Check the input */
1526 if(img==NULL) return STATUS_FAULT;
1527 if(img->status!=IMG_STATUS_OCCUPIED) return STATUS_FAULT;
1528 imgSetStatus(img, STATUS_FAULT);
1529 if(fname==NULL) return STATUS_FAULT;
1530 if(frame_index<0 || frame_index>img->dimt-1) return STATUS_FAULT;
1531 if(frame_to_read<1) return STATUS_FAULT;
1532
1533 /* Open the file */
1534 if((fp=fopen(fname, "rb")) == NULL) return STATUS_NOFILE;
1535
1536 /* Read main header */
1537 if((ret=ecat63ReadMainheader(fp, &main_header))) {
1538 fclose(fp); return STATUS_NOMAINHEADER;}
1539
1540 /* Read matrix list and nr and sort it */
1541 ecat63InitMatlist(&mlist);
1542 ret=ecat63ReadMatlist(fp, &mlist, IMG_TEST-1);
1543 if(ret) {fclose(fp); return STATUS_NOMATLIST;}
1544 if(mlist.matrixNr<=0) {
1545 fclose(fp); ecat63EmptyMatlist(&mlist); return STATUS_INVALIDMATLIST;}
1546 /* Make sure that plane and frame numbers are continuous */
1547 ecat63GatherMatlist(&mlist, 1, 1, 1, 1);
1549
1550 /* Calculate and check the size of one data matrix */
1551 ret=ecat63GetMatrixBlockSize(&mlist, &blkNr);
1552 if(ret) {ecat63EmptyMatlist(&mlist); fclose(fp); return ret;}
1553 /* And allocate memory for the raw data blocks */
1554 if(IMG_TEST>6) printf("allocating memory for %d blocks\n", blkNr);
1555 mdata=(char*)malloc((size_t)blkNr*MatBLKSIZE);
1556 if(mdata==NULL) {
1557 fclose(fp); ecat63EmptyMatlist(&mlist); return STATUS_NOMEMORY;}
1558
1559 /* Read all matrices that belong to the required frame */
1560 ret=0; seqplane=-1;
1561 for(int m=0; m<mlist.matrixNr; m++) if(mlist.matdir[m].matstat==1) {
1562 /* get frame and plane */
1563 mat_numdoc(mlist.matdir[m].matnum, &matval);
1564 plane=matval.plane;
1565 if(main_header.num_frames>=main_header.num_gates) frame=matval.frame;
1566 else frame=matval.gate; /* printf("frame=%d plane=%d\n", frame, plane); */
1567 if(frame!=frame_to_read) continue; /* do not process other frames */
1568 seqplane++;
1569 if(img->planeNumber[seqplane]<1) {
1570 img->planeNumber[seqplane]=plane;
1571 } else if(img->planeNumber[seqplane]!=plane) {
1572 fclose(fp); ecat63EmptyMatlist(&mlist); free(mdata);
1573 return STATUS_MISSINGMATRIX;
1574 }
1575
1576 /* Read the subheader */
1577 if(IMG_TEST>4) printf("reading subheader for matrix %d,%d\n", frame, plane);
1578 if(main_header.file_type==IMAGE_DATA) {
1579 ret=ecat63ReadImageheader(fp, mlist.matdir[m].strtblk, &image_header, IMG_TEST-10, NULL);
1580 scale=image_header.quant_scale;
1581 if(image_header.ecat_calibration_fctr>0.0 &&
1582 fabs(main_header.calibration_factor/image_header.ecat_calibration_fctr - 1.0)>0.0001)
1583 scale*=image_header.ecat_calibration_fctr;
1584 local_data_type=image_header.data_type;
1585 img->start[frame_index]=image_header.frame_start_time/1000.;
1586 img->end[frame_index]=img->start[frame_index]+image_header.frame_duration/1000.;
1587 img->mid[frame_index]=0.5*(img->start[frame_index]+img->end[frame_index]);
1588 if(image_header.decay_corr_fctr>1.0) {
1589 img->decayCorrFactor[frame_index]=image_header.decay_corr_fctr;
1591 } else {
1592 img->decayCorrFactor[frame_index]=0.0;
1594 }
1595 img->xform[0]=NIFTI_XFORM_UNKNOWN; // qform
1596 img->xform[1]=NIFTI_XFORM_SCANNER_ANAT; // sform
1597 img->quatern[6]=img->sizex; img->quatern[9]=img->sizex;
1598 img->quatern[11]=img->sizey; img->quatern[13]=img->sizey;
1599 img->quatern[16]=img->sizez; img->quatern[17]=img->sizez;
1600 } else if(main_header.file_type==RAW_DATA) {
1601 ret=ecat63ReadScanheader(fp, mlist.matdir[m].strtblk, &scan_header, IMG_TEST-10, NULL);
1602 scale=scan_header.scale_factor;
1603 if(scan_header.loss_correction_fctr>0.0) scale*=scan_header.loss_correction_fctr;
1604 local_data_type=scan_header.data_type;
1605 img->start[frame_index]=scan_header.frame_start_time/1000.;
1606 img->end[frame_index]=img->start[frame_index]+scan_header.frame_duration/1000.;
1607 img->mid[frame_index]=0.5*(img->start[frame_index]+img->end[frame_index]);
1608 img->sampleDistance=10.0*scan_header.sample_distance;
1609 img->prompts[frame_index]=(float)scan_header.prompts;
1610 img->randoms[frame_index]=(float)scan_header.delayed;
1611 } else if(main_header.file_type==ATTN_DATA) {
1612 ret=ecat63ReadAttnheader(fp, mlist.matdir[m].strtblk, &attn_header, IMG_TEST-10, NULL);
1613 img->sampleDistance=10.0*attn_header.sample_distance;
1614 scale=attn_header.scale_factor;
1615 local_data_type=attn_header.data_type;
1616 } else if(main_header.file_type==NORM_DATA) {
1617 ret=ecat63ReadNormheader(fp, mlist.matdir[m].strtblk, &norm_header, IMG_TEST-10, NULL);
1618 scale=norm_header.scale_factor;
1619 local_data_type=norm_header.data_type;
1620 } else
1621 img->_dataType=-1;
1622 if(ret) {
1623 ecat63EmptyMatlist(&mlist); free(mdata); fclose(fp);
1624 return STATUS_NOSUBHEADER;
1625 }
1626 if(main_header.calibration_factor>0.0) scale*=main_header.calibration_factor;
1627 if(IMG_TEST>6) printf("scale=%e datatype=%d\n", scale, img->_dataType);
1628 /* Read the pixel data */
1629 if(IMG_TEST>4) printf("reading matrix data\n");
1630 ret=ecat63ReadMatdata(fp, mlist.matdir[m].strtblk+1,
1631 mlist.matdir[m].endblk-mlist.matdir[m].strtblk,
1632 mdata, local_data_type);
1633 if(ret) {
1634 ecat63EmptyMatlist(&mlist); free(mdata); fclose(fp);
1635 return STATUS_MISSINGMATRIX;
1636 }
1637 if(IMG_TEST>4) printf("converting matrix data to floats\n");
1638 mptr=mdata;
1639 if(local_data_type==BYTE_TYPE) {
1640 for(int yi=0; yi<img->dimy; yi++)
1641 for(int xi=0; xi<img->dimx; xi++, mptr++) {
1642 img->m[seqplane][yi][xi][frame_index]=scale*(float)(*mptr);
1643 }
1644 } else if(local_data_type==VAX_I2 || local_data_type==SUN_I2) {
1645 for(int yi=0; yi<img->dimy; yi++)
1646 for(int xi=0; xi<img->dimx; xi++, mptr+=2) {
1647 sptr=(short int*)mptr;
1648 img->m[seqplane][yi][xi][frame_index]=scale*(float)(*sptr);
1649 }
1650 } else if(local_data_type==VAX_I4 || local_data_type==SUN_I4) {
1651 for(int yi=0; yi<img->dimy; yi++)
1652 for(int xi=0; xi<img->dimx; xi++, mptr+=4) {
1653 /*iptr=(int*)mptr;*/
1654 img->m[seqplane][yi][xi][frame_index]=1.0; /* scale*(float)(*iptr); */
1655 }
1656 } else if(local_data_type==VAX_R4 || local_data_type==IEEE_R4) {
1657 for(int yi=0; yi<img->dimy; yi++)
1658 for(int xi=0; xi<img->dimx; xi++, mptr+=4) {
1659 memcpy(&img->m[seqplane][yi][xi][frame_index], mptr, 4);
1660 img->m[seqplane][yi][xi][frame_index]*=scale;
1661 }
1662 }
1663 } /* next matrix */
1664 if(IMG_TEST>3) printf("end of matrices.\n");
1665
1666 free(mdata);
1667 ecat63EmptyMatlist(&mlist);
1668 fclose(fp);
1669
1670 /* seqplane is <0 only if this frame did not exist at all; return -1 */
1671 if(IMG_TEST>4) printf("last_seqplane := %d.\n", seqplane);
1672 if(seqplane<0) {imgSetStatus(img, STATUS_NOMATRIX); return STATUS_NOMATRIX;}
1673
1674 /* check that correct number of planes was read */
1675 if(seqplane+1 != img->dimz) {
1676 imgSetStatus(img, STATUS_MISSINGMATRIX);
1677 return STATUS_MISSINGMATRIX;
1678 }
1679
1680 /* For saving IMG data, only 2-byte VAX data types are allowed,
1681 so change it now */
1682 if(img->_dataType==VAX_I4 || img->_dataType==VAX_R4) img->_dataType=VAX_I2;
1683
1684 /* All went well */
1685 imgSetStatus(img, STATUS_OK);
1686 return STATUS_OK;
1687}
int ecat63GatherMatlist(MATRIXLIST *ml, short int do_planes, short int do_frames, short int do_gates, short int do_beds)
Definition ecat63ml.c:510
int ecat63GetMatrixBlockSize(MATRIXLIST *mlist, int *blk_nr)
Definition ecat63ml.c:366
void ecat63SortMatlistByFrame(MATRIXLIST *ml)
Definition ecat63ml.c:297

Referenced by imgReadEcat63FirstFrame(), and imgReadFrame().

◆ imgReadEcat63Header()

int imgReadEcat63Header ( const char * fname,
IMG * img )

Fill IMG structure header information from an image or sinogram file in ECAT 6.3 format. Information concerning separate frames or planes is not filled.

Note
ECAT 6.3 files do not have a magic number, therefore, do not use this function to determine if your file is in this format, at least test all other possible formats before calling this.
Parameters
fnameimage or sinogram filename
imgpointer to initialized IMG structure
Returns
errstatus, which is STATUS_OK (0) when call was successful, and >0 in case of an error.

Definition at line 1318 of file img_e63.c.

1321 {
1322 int m, blkNr=0, ret=0;
1323 int frameNr, planeNr;
1324 FILE *fp;
1325 ECAT63_mainheader main_header;
1326 MATRIXLIST mlist;
1327 ECAT63_imageheader image_header;
1328 ECAT63_scanheader scan_header;
1329 ECAT63_attnheader attn_header;
1330 ECAT63_normheader norm_header;
1331
1332
1333 if(IMG_TEST>0) printf("\nimgReadEcat63Header(%s, *img)\n", fname);
1334
1335 /* Check the arguments */
1336 if(img==NULL) return STATUS_FAULT;
1337 if(img->status!=IMG_STATUS_INITIALIZED) return STATUS_FAULT;
1338 imgSetStatus(img, STATUS_FAULT);
1339 if(fname==NULL) return STATUS_FAULT;
1340
1341 /* Open the file */
1342 if((fp=fopen(fname, "rb")) == NULL) return STATUS_NOFILE;
1343
1344 /* Read main header */
1345 if((ret=ecat63ReadMainheader(fp, &main_header))) {
1346 fclose(fp); return STATUS_NOMAINHEADER;}
1347 /* Check if file_type is supported */
1348 if(imgEcat63Supported(&main_header)==0) {
1349 fclose(fp); return STATUS_UNSUPPORTED;}
1350 /* Copy main header information into IMG; sets also img.type */
1351 imgGetEcat63MHeader(img, &main_header);
1352 if(IMG_TEST>7) printf("img.type := %d\n", img->type);
1353 img->_fileFormat=imgGetEcat63Fileformat(&main_header);
1354 if(IMG_TEST>7) printf("img._fileFormat := %d\n", img->_fileFormat);
1355 if(img->_fileFormat==IMG_UNKNOWN) {fclose(fp); return STATUS_UNSUPPORTED;}
1356
1357 /* Read matrix list and nr and sort it */
1358 ecat63InitMatlist(&mlist);
1359 ret=ecat63ReadMatlist(fp, &mlist, IMG_TEST-1);
1360 if(ret) {fclose(fp); return STATUS_NOMATLIST;}
1361 if(mlist.matrixNr<=0) {
1362 ecat63EmptyMatlist(&mlist); fclose(fp); return STATUS_INVALIDMATLIST;}
1363 /* Make sure that plane and frame numbers are continuous */
1364 ecat63GatherMatlist(&mlist, 1, 1, 1, 1);
1365 ecat63SortMatlistByPlane(&mlist); // otherwise frameNr can't be computed below
1366 /* Trim extra frames */
1367 if(main_header.num_frames>0)
1368 /*del_nr=*/ecat63DeleteLateFrames(&mlist, main_header.num_frames);
1369 /* Get plane and frame numbers and check that volume is full */
1370 ret=ecat63GetPlaneAndFrameNr(&mlist, &main_header, &planeNr, &frameNr);
1371 if(ret) {ecat63EmptyMatlist(&mlist); fclose(fp); return ret;}
1372 img->dimz=planeNr;
1373 img->dimt=frameNr;
1374 /* Calculate and check the size of one data matrix */
1375 ret=ecat63GetMatrixBlockSize(&mlist, &blkNr);
1376 if(ret) {ecat63EmptyMatlist(&mlist); fclose(fp); return ret;}
1377
1378 /* Read one subheader */
1379 if(IMG_TEST>5) printf("main_header.file_type := %d\n", main_header.file_type);
1380 m=0;
1381 switch(main_header.file_type) {
1382 case IMAGE_DATA:
1383 ret=ecat63ReadImageheader(fp, mlist.matdir[m].strtblk, &image_header, IMG_TEST-10, NULL);
1384 break;
1385 case RAW_DATA:
1386 ret=ecat63ReadScanheader(fp, mlist.matdir[m].strtblk, &scan_header, IMG_TEST-10, NULL);
1387 break;
1388 case ATTN_DATA:
1389 ret=ecat63ReadAttnheader(fp, mlist.matdir[m].strtblk, &attn_header, IMG_TEST-10, NULL);
1390 break;
1391 case NORM_DATA:
1392 ret=ecat63ReadNormheader(fp, mlist.matdir[m].strtblk, &norm_header, IMG_TEST-10, NULL);
1393 break;
1394 default: ret=-1;
1395 }
1396 /* Free locally allocated memory and close the file */
1397 ecat63EmptyMatlist(&mlist); fclose(fp);
1398 /* Check whether subheader was read */
1399 if(ret) return STATUS_NOSUBHEADER;
1400
1401 /* Get the following information in the subheader:
1402 dimensions x, y and z; datatype;
1403 image decay correction on/off, zoom, and pixel size;
1404 sinogram sample distance.
1405 */
1406 switch(main_header.file_type) {
1407 case IMAGE_DATA:
1408 img->dimx=image_header.dimension_1; img->dimy=image_header.dimension_2;
1409 if(img->unit<1) imgUnitFromEcat(img, image_header.quant_units);
1410 img->_dataType=image_header.data_type;
1411 img->zoom=image_header.recon_scale;
1412 if(image_header.decay_corr_fctr>1.0) img->decayCorrection=IMG_DC_CORRECTED;
1414 img->sizex=img->sizey=10.*image_header.pixel_size;
1415 if(img->sizez<10.*image_header.slice_width)
1416 img->sizez=10.*image_header.slice_width;
1417 break;
1418 case RAW_DATA:
1419 img->dimx=scan_header.dimension_1; img->dimy=scan_header.dimension_2;
1420 img->type=IMG_TYPE_RAW;
1421 img->_dataType=scan_header.data_type;
1423 img->sampleDistance=10.0*scan_header.sample_distance;
1424 break;
1425 case ATTN_DATA:
1426 img->dimx=attn_header.dimension_1; img->dimy=attn_header.dimension_2;
1427 img->type=IMG_TYPE_ATTN;
1429 img->_dataType=attn_header.data_type;
1430 break;
1431 case NORM_DATA:
1432 img->dimx=norm_header.dimension_1; img->dimy=norm_header.dimension_2;
1433 img->type=IMG_TYPE_RAW;
1435 img->_dataType=norm_header.data_type;
1436 break;
1437 }
1438
1439 /* For saving IMG data, only 2-byte VAX data types are allowed,
1440 so change it now */
1441 if(img->_dataType==VAX_I4 || img->_dataType==VAX_R4) img->_dataType=VAX_I2;
1442
1443 imgSetStatus(img, STATUS_OK);
1444 return STATUS_OK;
1445}
int ecat63GetPlaneAndFrameNr(MATRIXLIST *mlist, ECAT63_mainheader *h, int *plane_nr, int *frame_nr)
Definition ecat63ml.c:400
void ecat63SortMatlistByPlane(MATRIXLIST *ml)
Definition ecat63ml.c:271
int imgGetEcat63Fileformat(ECAT63_mainheader *h)
Definition img_e63.c:1288
int imgEcat63Supported(ECAT63_mainheader *h)
Definition img_e63.c:1155
void imgGetEcat63MHeader(IMG *img, ECAT63_mainheader *h)
Definition img_e63.c:1172

Referenced by imgFormatDetermine(), imgReadEcat63FirstFrame(), imgReadHeader(), and imgWriteEcat63Frame().

◆ imgSetEcat63MHeader()

void imgSetEcat63MHeader ( IMG * img,
ECAT63_mainheader * h )

Copy information from IMG struct into ECAT 6.3 main header

Parameters
imgsource image structure
htarget Ecat 6.3 main header

Definition at line 1223 of file img_e63.c.

1224{
1225 struct tm tm;
1226 memset(&tm, 0, sizeof(struct tm));
1227
1228 if(IMG_TEST>0) printf("imgSetEcat63MHeader()\n");
1229 if(IMG_TEST>2) {
1230 char buf[32];
1231 if(!ctime_r_int(&img->scanStart, buf)) strcpy(buf, "1900-01-01 00:00:00");
1232 fprintf(stdout, " scan_start_time := %s\n", buf);
1233 }
1234 h->sw_version=2;
1235 h->num_planes=img->dimz;
1236 h->num_frames=img->dimt;
1237 h->num_gates=1;
1238 h->num_bed_pos=1;
1240 else h->file_type=RAW_DATA;
1241 h->data_type=VAX_I2;
1242 if(img->scanner>0) h->system_type=img->scanner;
1244 h->calibration_factor=1.0;
1245 h->axial_fov=img->axialFOV/10.0;
1246 h->transaxial_fov=img->transaxialFOV/10.0;
1247 h->plane_separation=img->sizez/10.0;
1249 strncpy(h->radiopharmaceutical, img->radiopharmaceutical, 32);
1250 if(gmtime_r(&img->scanStart, &tm)!=NULL) {
1251 h->scan_start_year=tm.tm_year+1900;
1252 h->scan_start_month=tm.tm_mon+1;
1253 h->scan_start_day=tm.tm_mday;
1254 h->scan_start_hour=tm.tm_hour;
1255 h->scan_start_minute=tm.tm_min;
1256 h->scan_start_second=tm.tm_sec;
1257 if(IMG_TEST>2) {
1258 printf(" img->scanStart := %ld\n", (long int)img->scanStart);
1259 printf(" -> tm_year := %d\n", tm.tm_year);
1260 printf(" -> tm_hour := %d\n", tm.tm_hour);
1261 }
1262 } else {
1263 h->scan_start_year=1900;
1264 h->scan_start_month=1;
1265 h->scan_start_day=1;
1266 h->scan_start_hour=0;
1267 h->scan_start_minute=0;
1268 h->scan_start_second=0;
1269 if(IMG_TEST>0) printf("invalid scan_start_time in IMG\n");
1270 }
1272 strcpy(h->isotope_code, imgIsotope(img));
1273 strlcpy(h->study_name, img->studyNr, 12);
1274 strcpy(h->patient_name, img->patientName);
1275 strcpy(h->patient_id, img->patientID);
1278}

Referenced by imgWriteEcat63Frame().

◆ imgSetEcat63SHeader()

void imgSetEcat63SHeader ( IMG * img,
void * h )

Copies Ecat6.3 image or scan sub header information from IMG struct.

Parameters
imgsource image structure
htarget sub header structure

Definition at line 1872 of file img_e63.c.

1872 {
1873 ECAT63_imageheader *image_header;
1874 ECAT63_scanheader *scan_header;
1875
1876 if(img->type==IMG_TYPE_RAW) {
1877 scan_header=(ECAT63_scanheader*)h;
1878 scan_header->data_type=VAX_I2;
1879 scan_header->dimension_1=img->dimx;
1880 scan_header->dimension_2=img->dimy;
1881 scan_header->frame_duration_sec=1;
1882 scan_header->scale_factor=1.0;
1883 scan_header->frame_start_time=0;
1884 scan_header->frame_duration=1000;
1885 scan_header->loss_correction_fctr=1.0;
1886 scan_header->sample_distance=img->sampleDistance/10.0;
1887 } else {
1888 image_header=(ECAT63_imageheader*)h;
1889 image_header->data_type=VAX_I2;
1890 image_header->num_dimensions=2;
1891 image_header->dimension_1=img->dimx;
1892 image_header->dimension_2=img->dimy;
1893 image_header->recon_scale=img->zoom;
1894 image_header->quant_scale=1.0;
1895 image_header->slice_width=img->sizez/10.;
1896 image_header->pixel_size=img->sizex/10.;
1897 image_header->frame_start_time=0;
1898 image_header->frame_duration=1000;
1899 image_header->plane_eff_corr_fctr=1.0; // now we don't know which plane
1900 image_header->decay_corr_fctr=1.0;
1901 image_header->loss_corr_fctr=1.0;
1902 image_header->ecat_calibration_fctr=1.0;
1903 image_header->well_counter_cal_fctr=1.0;
1904 image_header->quant_units=imgUnitToEcat6(img);
1905 }
1906}

Referenced by imgWriteEcat63Frame().

◆ imgWriteEcat63Frame()

int imgWriteEcat63Frame ( const char * fname,
int frame_to_write,
IMG * img,
int frame_index )

Write one PET frame from IMG data struct into ECAT 6.3 image or sinogram file; format is specified in IMG struct. This function can be called repeatedly to write all frames one at a time to conserve memory. However, file with just mainheader and matrix list without any previous frame is not accepted.

Parameters
fnamename of file where IMG contents will be written. If file does not exist, it is created. Make sure to delete existing file, unless you want to add data
frame_to_writePET frame number (1..frameNr) which will be written: If set to 0, frame data will be written to an existing or new PET file as a new frame, never overwriting existing data. If >0, then frame data is written as specified frame number, overwriting any data existing with the same frame number
imgpointer to the IMG data struct
frame_indexIMG frame index (0..dimt-1) which will be written
Returns
errstatus, which is STATUS_OK (0) when call was successful, and >0 in case of an error.

Definition at line 1711 of file img_e63.c.

1713 {
1714 IMG test_img;
1715 int ret=0, matrixId;
1716 ECAT63_mainheader main_header;
1717 ECAT63_imageheader image_header;
1718 ECAT63_scanheader scan_header;
1719 void *sub_header=NULL;
1720 FILE *fp;
1721 float *fdata=NULL, *fptr;
1722
1723
1724 if(IMG_TEST>0) printf("\nimgWriteEcat63Frame(%s, %d, *img, %d)\n",
1725 fname, frame_to_write, frame_index);
1727 if(IMG_TEST>4) {
1728 char buf[32];
1729 if(!ctime_r_int(&img->scanStart, buf)) strcpy(buf, "1900-01-01 00:00:00");
1730 fprintf(stdout, " scan_start_time := %s\n", buf);
1731 }
1732
1733 /*
1734 * Check the input
1735 */
1736 if(fname==NULL) return STATUS_FAULT;
1737 if(img==NULL) return STATUS_FAULT;
1738 if(img->status!=IMG_STATUS_OCCUPIED) return STATUS_FAULT;
1739 if(frame_to_write<0) return STATUS_FAULT;
1740 if(frame_index<0 || frame_index>=img->dimt) return STATUS_FAULT;
1741 if(img->_fileFormat!=IMG_E63) return STATUS_FAULT;
1742
1743 /* Check image size */
1744 if(img->dimt>511 || img->dimz>255 || img->dimx>1024 || img->dimy>1024) {
1745 strcpy(ecat63errmsg, "too large matrix size"); return(1);}
1746
1747 /* Initiate headers */
1748 memset(&main_header, 0, sizeof(ECAT63_mainheader));
1749 memset(&image_header,0, sizeof(ECAT63_imageheader));
1750 memset(&scan_header, 0, sizeof(ECAT63_scanheader));
1751 imgInit(&test_img);
1752
1753
1754 /*
1755 * If file does not exist, then create it with new mainheader,
1756 * and if it does exist, then read and check header information.
1757 * Create or edit mainheader to contain correct frame nr.
1758 * In any case, leave file pointer open for write.
1759 */
1760 if(access(fname, 0) == -1) { /* file does not exist*/
1761
1762 if(IMG_TEST>1) printf(" new file\n");
1763 /* Set main header */
1764 imgSetEcat63MHeader(img, &main_header);
1765 if(IMG_TEST>6) {
1766 ecat63PrintMainheader(&main_header, stdout);
1767 }
1768 if(frame_to_write==0) frame_to_write=1;
1769 main_header.num_frames=frame_to_write;
1770
1771 /* Open file, write main header and initiate matrix list */
1772 fp=ecat63Create(fname, &main_header); if(fp==NULL) return STATUS_NOWRITEPERM;
1773
1774 } else { /* file exists*/
1775
1776 if(IMG_TEST>1) printf(" existing file\n");
1777 /* Read header information for checking */
1778 ret=imgReadEcat63Header(fname, &test_img); if(ret!=0) return ret;
1779 if(IMG_TEST>1) {
1780 char buf[32];
1781 if(!ctime_r_int(&test_img.scanStart, buf))
1782 strcpy(buf, "1900-01-01 00:00:00");
1783 fprintf(stdout, "scan_start_time := %s\n", buf);
1784 }
1785 /* Check that file format is the same */
1786 if(img->_fileFormat!=test_img._fileFormat || img->type!=test_img.type)
1787 return STATUS_WRONGFILETYPE;
1788 /* Check that matrix sizes are the same */
1789 if(img->dimz!=test_img.dimz || img->dimx!=test_img.dimx ||
1790 img->dimy!=test_img.dimy)
1791 return STATUS_VARMATSIZE;
1792 imgEmpty(&test_img);
1793
1794 /* Open the file for read and write */
1795 if((fp=fopen(fname, "r+b"))==NULL) return STATUS_NOWRITEPERM;
1796
1797 /* Read the mainheader, set new frame number, and write it back */
1798 if((ret=ecat63ReadMainheader(fp, &main_header))!=0)
1799 return STATUS_NOMAINHEADER;
1800 if(frame_to_write==0) frame_to_write=main_header.num_frames+1;
1801 if(main_header.num_frames<frame_to_write)
1802 main_header.num_frames=frame_to_write;
1803 if((ret=ecat63WriteMainheader(fp, &main_header))!=0)
1804 return STATUS_NOWRITEPERM;
1805 if(IMG_TEST>0) {
1806 char tmp[32];
1807 printf(" scan_start_time := %s\n", ecat63ScanstarttimeInt(&main_header, tmp));
1808 }
1809
1810 }
1811 if(IMG_TEST>2) {
1812 printf("frame_to_write := %d\n", frame_to_write);
1813 }
1814
1815 /* Allocate memory for matrix float data */
1816 size_t pxlNr=img->dimx*img->dimy*img->dimz;
1817 fdata=(float*)malloc(pxlNr*sizeof(float));
1818 if(fdata==NULL) {fclose(fp); return STATUS_NOMEMORY;}
1819
1820 /* Set matrix subheader */
1821 if(img->type==IMG_TYPE_RAW) sub_header=(void*)&scan_header;
1822 else if(img->type==IMG_TYPE_IMAGE) sub_header=&image_header;
1823 else {fclose(fp); return STATUS_FAULT;}
1824 imgSetEcat63SHeader(img, sub_header);
1825
1826 /* Copy matrix pixel values to fdata */
1827 fptr=fdata;
1828 for(int zi=0; zi<img->dimz; zi++)
1829 for(int yi=0; yi<img->dimy; yi++)
1830 for(int xi=0; xi<img->dimx; xi++)
1831 *fptr++=img->m[zi][yi][xi][frame_index];
1832
1833 /* Write subheader and data, and set the rest of subheader contents */
1834 fptr=fdata; ret=-1;
1835 for(int zi=0; zi<img->dimz; zi++, fptr+=img->dimx*img->dimy) {
1836 /* Create new matrix id (i.e. matnum) */
1837 matrixId=mat_numcod(frame_to_write, img->planeNumber[zi], 1, 0, 0);
1838 if(img->type==IMG_TYPE_RAW) {
1839 scan_header.frame_start_time=
1840 (int)temp_roundf(1000.*img->start[frame_index]);
1841 scan_header.frame_duration=
1842 (int)temp_roundf(1000.*(img->end[frame_index]-img->start[frame_index]));
1843 scan_header.prompts=temp_roundf(img->prompts[frame_index]);
1844 scan_header.delayed=temp_roundf(img->randoms[frame_index]);
1845 ret=ecat63WriteScanMatrix(fp, matrixId, &scan_header, fptr);
1846 } else {
1847 image_header.frame_start_time=
1848 (int)temp_roundf(1000.*img->start[frame_index]);
1849 image_header.frame_duration=
1850 (int)temp_roundf(1000.*(img->end[frame_index]-img->start[frame_index]));
1852 image_header.decay_corr_fctr=img->decayCorrFactor[frame_index];
1853 else
1854 image_header.decay_corr_fctr=0.0;
1855 ret=ecat63WriteImageMatrix(fp, matrixId, &image_header, fptr);
1856 }
1857 if(ret!=0) break;
1858 } /* next plane*/
1859 free(fdata); fclose(fp);
1860 if(ret) return STATUS_DISKFULL;
1861
1862 return STATUS_OK;
1863}
int ecat63WriteScanMatrix(FILE *fp, int matnum, ECAT63_scanheader *h, float *fdata)
Definition ecat63w.c:781
int ecat63WriteImageMatrix(FILE *fp, int matnum, ECAT63_imageheader *h, float *fdata)
Definition ecat63w.c:697
void imgInit(IMG *image)
Definition img.c:60
void imgSetEcat63SHeader(IMG *img, void *h)
Definition img_e63.c:1872
void imgSetEcat63MHeader(IMG *img, ECAT63_mainheader *h)
Definition img_e63.c:1223

Referenced by imgWriteFrame().