00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048 #include <stdio.h>
00049 #include <stdlib.h>
00050 #include <unistd.h>
00051 #include <math.h>
00052 #include <string.h>
00053 #include <time.h>
00054
00055 #include "petc99.h"
00056 #include "swap.h"
00057 #include "halflife.h"
00058
00059 #include "include/img.h"
00060 #include "include/analyze.h"
00061 #include "include/imgmax.h"
00062 #include "include/imgdecay.h"
00063 #include "include/sif.h"
00064 #include "include/imgfile.h"
00065
00066
00067
00083 int imgReadAnalyze(const char *dbname, IMG *img) {
00084 FILE *fp;
00085 int ret, fi, pi, xi, yi;
00086 float *fdata=NULL, *fptr;
00087 ANALYZE_DSR dsr;
00088 char datfile[FILENAME_MAX], hdrfile[FILENAME_MAX], siffile[FILENAME_MAX];
00089 char buf[128], *cptr;
00090 int dimNr, dimx, dimy, dimz=1, dimt=1, pxlNr=0;
00091 SIF sif;
00092 struct tm *st;
00093
00094
00095 if(IMG_TEST) printf("imgReadAnalyze(%s, *img)\n", dbname);
00096
00097
00098 imgSetStatus(img, STATUS_OK);
00099 if(img==NULL || img->status!=IMG_STATUS_INITIALIZED) {
00100 imgSetStatus(img, STATUS_FAULT); return(2);}
00101 if(dbname==NULL || !dbname[0]) {imgSetStatus(img, STATUS_FAULT); return(1);}
00102
00103
00104 if(anaExists(dbname)==0) {
00105
00106 strcpy(datfile, dbname); cptr=strrchr(datfile, '.');
00107 if(cptr!=NULL && (strcmp(cptr, ".img")==0 || strcmp(cptr, ".hdr")==0)) {
00108 *cptr=(char)0; strcpy(hdrfile, datfile);
00109 if(anaExists(datfile)==0) {
00110 imgSetStatus(img, STATUS_NOHEADERFILE); return(3);}
00111 strcat(datfile, ".img"); strcat(hdrfile, ".hdr");
00112 } else {
00113 imgSetStatus(img, STATUS_NOHEADERFILE); return(3);
00114 }
00115 } else {
00116
00117 strcpy(datfile, dbname); strcat(datfile, ".img");
00118 strcpy(hdrfile, dbname); strcat(hdrfile, ".hdr");
00119 }
00120
00121
00122 ret=anaReadHeader(hdrfile, &dsr);
00123 if(ret) {
00124 if(ret==1) imgSetStatus(img, STATUS_FAULT);
00125 else if(ret==2) imgSetStatus(img, STATUS_NOHEADERFILE);
00126 else imgSetStatus(img, STATUS_UNSUPPORTED);
00127 return(3);
00128 }
00129 if(IMG_TEST) anaPrintHeader(&dsr, stdout);
00130
00131
00132 if(IMG_TEST) fprintf(stdout, "reading image data %s\n", datfile);
00133 if((fp=fopen(datfile, "rb")) == NULL) {imgSetStatus(img, STATUS_NOIMGDATA); return(5);}
00134
00135
00136
00137 dimNr=dsr.dime.dim[0];
00138 if(dimNr<2) {fclose(fp); imgSetStatus(img, STATUS_INVALIDHEADER); return(4);}
00139 dimx=dsr.dime.dim[1]; dimy=dsr.dime.dim[2];
00140 if(dimNr>2) {dimz=dsr.dime.dim[3]; if(dimNr>3) dimt=dsr.dime.dim[4];}
00141 pxlNr=dimx*dimy*dimz;
00142 if(pxlNr<1) {fclose(fp); imgSetStatus(img, STATUS_INVALIDHEADER); return(4);}
00143
00144 ret=imgAllocate(img, dimz, dimy, dimx, dimt);
00145 if(ret) {fclose(fp); imgSetStatus(img, STATUS_NOMEMORY); return(11);}
00146
00147 img->type=IMG_TYPE_IMAGE;
00148 strncpy(img->studyNr, dsr.hist.patient_id, MAX_STUDYNR_LEN);
00149 strcpy(img->patientName, dsr.hist.patient_id);
00150 img->sizex=dsr.dime.pixdim[1];
00151 img->sizey=dsr.dime.pixdim[2];
00152 img->sizez=dsr.dime.pixdim[3];
00153
00154 if(dsr.dime.funused3>1.E-5) img->isotopeHalflife=dsr.dime.funused3;
00155 for(pi=0; pi<dimz; pi++) img->planeNumber[pi]=pi+1;
00156 if(dsr.little) img->_fileFormat=IMG_ANA_L; else img->_fileFormat=IMG_ANA;
00157
00158 if(strstr(dsr.hist.descrip, "Decay corrected.")!=NULL)
00159 img->decayCorrected=1;
00160 else if(strstr(dsr.hist.descrip, "No decay correction.")!=NULL)
00161 img->decayCorrected=0;
00162 else
00163 img->decayCorrected=1;
00164
00165
00166 fdata=malloc(pxlNr*sizeof(float));
00167 if(fdata==NULL) {fclose(fp); imgSetStatus(img, STATUS_NOMEMORY); return(12);}
00168
00169
00170 for(fi=0; fi<dimt; fi++) {
00171 fptr=fdata;
00172 ret=anaReadImagedata(fp, &dsr, fi+1, fptr);
00173 if(ret) {free(fdata); fclose(fp); imgSetStatus(img, STATUS_NOIMGDATA); return(7);}
00174
00175 fptr=fdata;
00176 if(anaFlipping()==0) {
00177 for(pi=0; pi<img->dimz; pi++)
00178 for(yi=dimy-1; yi>=0; yi--)
00179 for(xi=dimx-1; xi>=0; xi--)
00180 img->m[pi][yi][xi][fi]=*fptr++;
00181 } else {
00182 for(pi=dimz-1; pi>=0; pi--)
00183 for(yi=dimy-1; yi>=0; yi--)
00184 for(xi=dimx-1; xi>=0; xi--)
00185 img->m[pi][yi][xi][fi]=*fptr++;
00186 }
00187 }
00188 free(fdata);
00189 fclose(fp);
00190
00191
00192
00193 strcpy(siffile, dbname); strcat(siffile, ".sif");
00194 if(access(siffile, 0) == -1) {
00195 strcpy(siffile, datfile); strcat(siffile, ".sif");
00196 }
00197
00198 if(IMG_TEST) printf("reading SIF file %s\n", siffile);
00199 if(access(siffile, 0) == -1) {
00200 if(IMG_TEST) printf(" No SIF file; therefore unknown frame times.\n");
00201 return(0);
00202 }
00203
00204 sifInit(&sif); ret=sifRead(siffile, &sif);
00205 if(ret) {imgSetStatus(img, STATUS_NOSIFDATA); return(21);}
00206
00207 if(sif.scantime>0) {
00208 st=localtime(&sif.scantime);
00209 if(st!=NULL) strftime(buf, 128, "%Y-%m-%d %H:%M:%S", st);
00210 else strcpy(buf, "1900-01-01 00:00:00");
00211 img->scanStart=sif.scantime;
00212 }
00213
00214 if(sif.frameNr!=img->dimt) {
00215 imgSetStatus(img, STATUS_WRONGSIFDATA); sifEmpty(&sif); return(22);}
00216 for(fi=0; fi<sif.frameNr; fi++) {
00217 img->start[fi]=sif.x1[fi]; img->end[fi]=sif.x2[fi];
00218 img->mid[fi]=0.5*(img->start[fi]+img->end[fi]);
00219 }
00220
00221 for(fi=0; fi<sif.frameNr; fi++) {
00222 img->prompts[fi]=sif.prompts[fi];
00223 img->randoms[fi]=sif.randoms[fi];
00224 }
00225
00226
00227 if(strlen(sif.isotope_name)>1) {
00228 img->isotopeHalflife=60.0*hlFromIsotope(sif.isotope_name);
00229 }
00230 sifEmpty(&sif);
00231
00232 return(0);
00233 }
00234
00235
00236
00253 int imgWriteAnalyze(const char *dbname, IMG *img) {
00254 FILE *fp;
00255 int ret, fi, pi, xi, yi, little;
00256 float g;
00257 ANALYZE_DSR dsr;
00258 char datfile[FILENAME_MAX], hdrfile[FILENAME_MAX];
00259 const char *cptr;
00260 int pxlNr=0;
00261 struct tm *st;
00262 short int *sdata, *sptr, smin, smax;
00263
00264
00265 if(IMG_TEST) printf("imgWriteAnalyze(%s, *img)\n", dbname);
00266
00267
00268 imgSetStatus(img, STATUS_OK);
00269 if(img==NULL || img->status!=IMG_STATUS_OCCUPIED) {
00270 imgSetStatus(img, STATUS_FAULT); return(2);}
00271 if(dbname==NULL || !dbname[0]) {imgSetStatus(img, STATUS_FAULT); return(1);}
00272
00273
00274 strcpy(datfile, dbname); strcat(datfile, ".img");
00275 strcpy(hdrfile, dbname); strcat(hdrfile, ".hdr");
00276
00277
00278
00279
00280
00281 if(img->_fileFormat==IMG_ANA_L) dsr.little=1; else dsr.little=0;
00282
00283 memset(&dsr.hk, 0, sizeof(ANALYZE_HEADER_KEY));
00284 memset(&dsr.dime, 0, sizeof(ANALYZE_HEADER_IMGDIM));
00285 memset(&dsr.hist, 0, sizeof(ANALYZE_HEADER_HISTORY));
00286 dsr.hk.sizeof_hdr=348;
00287 strcpy(dsr.hk.data_type, "");
00288 cptr=strrchr(dbname, '/'); if(cptr==NULL) cptr=strrchr(dbname, '\\');
00289 if(cptr!=NULL) cptr++; if(cptr==NULL) cptr=dbname;
00290 strncpy(dsr.hk.db_name, cptr, 17);
00291 dsr.hk.extents=16384;
00292 dsr.hk.regular='r';
00293
00294 dsr.dime.dim[0]=4;
00295 dsr.dime.dim[1]=img->dimx;
00296 dsr.dime.dim[2]=img->dimy;
00297 dsr.dime.dim[3]=img->dimz;
00298 dsr.dime.dim[4]=img->dimt;
00299 dsr.dime.datatype=ANALYZE_DT_SIGNED_SHORT;
00300 dsr.dime.bitpix=16;
00301 dsr.dime.pixdim[0]=0.0;
00302 dsr.dime.pixdim[1]=img->sizex;
00303 dsr.dime.pixdim[2]=img->sizey;
00304 dsr.dime.pixdim[3]=img->sizez;
00305 dsr.dime.pixdim[4]=0.0;
00306 dsr.dime.funused1=0.0;
00307
00308 dsr.dime.funused3=img->isotopeHalflife;
00309
00310 if(img->decayCorrected==1) strcpy(dsr.hist.descrip, "Decay corrected.");
00311 else strcpy(dsr.hist.descrip, "No decay correction.");
00312 strncpy(dsr.hist.scannum, img->studyNr, 10);
00313 st=localtime(&img->scanStart);
00314 if(st!=NULL) {
00315 strftime(dsr.hist.exp_date, 10, "%Y-%m-%d", st);
00316 strftime(dsr.hist.exp_time, 10, "%H:%M:%S", st);
00317 } else {
00318 strncpy(dsr.hist.exp_date, "1900-01-01", 10);
00319 strncpy(dsr.hist.exp_time, "00:00:00", 10);
00320 }
00321
00322
00323
00324
00325
00326 if(IMG_TEST) printf("scaling data to short ints\n");
00327 ret=imgMinMax(img, &dsr.dime.cal_min, &dsr.dime.cal_max);
00328 if(ret) {imgSetStatus(img, STATUS_FAULT); return(3);}
00329 if(fabs(dsr.dime.cal_min)>fabs(dsr.dime.cal_max)) g=fabs(dsr.dime.cal_min);
00330 else g=fabs(dsr.dime.cal_max);
00331 if(g<1E-20) g=1.0; else g=32767./g; dsr.dime.funused1=1.0/g;
00332 if(IMG_TEST) printf("min=%g max=%g scale_factor=%g\n",
00333 dsr.dime.cal_min, dsr.dime.cal_max, dsr.dime.funused1);
00334
00335
00336 pxlNr=(img->dimx)*(img->dimy)*(img->dimz);
00337 sdata=malloc(pxlNr*sizeof(short int));
00338 if(sdata==NULL) {
00339 imgSetStatus(img, STATUS_NOMEMORY);
00340 return 12;
00341 }
00342
00343
00344 if((fp=fopen(datfile, "wb")) == NULL) {
00345 imgSetStatus(img, STATUS_CANTWRITEIMGFILE);
00346 free(sdata);
00347 return 14;
00348 }
00349
00350
00351
00352 smin=smax=temp_roundf(g*img->m[0][0][0][0]);
00353 for(fi=0; fi<img->dimt; fi++) {
00354 sptr=sdata;
00355 if(anaFlipping()==0) {
00356 for(pi=0; pi<img->dimz; pi++)
00357 for(yi=img->dimy-1; yi>=0; yi--)
00358 for(xi=img->dimx-1; xi>=0; xi--) {
00359 *sptr=temp_roundf(g*img->m[pi][yi][xi][fi]);
00360 if(*sptr>smax) smax=*sptr; else if(*sptr<smin) smin=*sptr;
00361 sptr++;
00362 }
00363 } else {
00364 for(pi=img->dimz-1; pi>=0; pi--)
00365 for(yi=img->dimy-1; yi>=0; yi--)
00366 for(xi=img->dimx-1; xi>=0; xi--) {
00367 *sptr=temp_roundf(g*img->m[pi][yi][xi][fi]);
00368 if(*sptr>smax) smax=*sptr; else if(*sptr<smin) smin=*sptr;
00369 sptr++;
00370 }
00371 }
00372
00373 little=little_endian();
00374 if(little!=dsr.little)
00375 swabip(sdata, pxlNr*sizeof(short int));
00376
00377 if(fwrite(sdata, 2, pxlNr, fp) != pxlNr) {
00378 imgSetStatus(img, STATUS_CANTWRITEIMGFILE);
00379 free(sdata); fclose(fp);
00380 return 15;
00381 }
00382 }
00383
00384 fclose(fp);
00385 free(sdata);
00386
00387 if(IMG_TEST) printf("smin=%d smax=%d\n", smin, smax);
00388
00389
00390 dsr.dime.glmin=smin; dsr.dime.glmax=smax;
00391
00392
00393 ret=anaWriteHeader(hdrfile, &dsr);
00394 if(ret) {
00395 imgSetStatus(img, STATUS_CANTWRITEHEADERFILE);
00396 return 21;
00397 }
00398 imgSetStatus(img, STATUS_OK);
00399 return 0;
00400 }
00401
00402
00403
00414 int imgReadAnalyzeHeader(const char *dbname, IMG *img) {
00415 char hdrfile[FILENAME_MAX], siffile[FILENAME_MAX];
00416 ANALYZE_DSR ana_header;
00417 SIF sif;
00418 double f;
00419 int ret;
00420
00421 if(IMG_TEST) printf("\nimgReadAnalyzeHeader(%s, *img)\n", dbname);
00422
00423
00424 if(img==NULL) return STATUS_FAULT;
00425 if(img->status!=IMG_STATUS_INITIALIZED) return STATUS_FAULT;
00426 imgSetStatus(img, STATUS_FAULT);
00427 if(dbname==NULL) return STATUS_FAULT;
00428
00429
00430 ret=anaDatabaseExists(dbname, hdrfile, NULL, siffile);
00431 if(ret==0) return STATUS_NOFILE;
00432
00433
00434 ret=anaReadHeader(hdrfile, &ana_header);
00435 if(ret!=0) {
00436 if(IMG_TEST>1) printf("anaReadHeader() return value := %d\n", ret);
00437 if(ret==1) return STATUS_FAULT;
00438 else if(ret==2) return STATUS_NOHEADERFILE;
00439 else return STATUS_UNSUPPORTED;
00440 return(STATUS_FAULT);
00441 }
00442
00443 ret=imgGetAnalyzeHeader(img, &ana_header);
00444 if(ret!=0) {
00445 imgSetStatus(img, ret);
00446 return(ret);
00447 }
00448
00449
00450 if(!siffile[0]) {
00451 imgSetStatus(img, STATUS_OK);
00452 return STATUS_OK;
00453 }
00454
00455
00456 sifInit(&sif); ret=0;
00457 if(sifRead(siffile, &sif)!=0) return STATUS_OK;
00458
00459 img->scanStart=sif.scantime;
00460
00461 if(!img->studyNr[0] && strlen(sif.studynr)>1 )
00462 strncpy(img->studyNr, sif.studynr, MAX_STUDYNR_LEN);
00463
00464 f=hlFromIsotope(sif.isotope_name);
00465 if(img->isotopeHalflife<=0.0 && f>0.0) img->isotopeHalflife=60.0*f;
00466 sifEmpty(&sif);
00467
00468 return STATUS_OK;
00469 }
00470
00471
00472
00481 int imgGetAnalyzeHeader(IMG *img, ANALYZE_DSR *h) {
00482 int dimNr, dimx, dimy, dimz=1, dimt=1, pxlNr=0;
00483
00484 if(IMG_TEST) printf("\nimgGetAnalyzeHeader(*img, *dsr)\n");
00485
00486
00487
00488 if(img==NULL) return STATUS_FAULT;
00489 if(img->status!=IMG_STATUS_INITIALIZED && img->status!=IMG_STATUS_OCCUPIED)
00490 return STATUS_FAULT;
00491 imgSetStatus(img, STATUS_FAULT);
00492 if(h==NULL) return STATUS_FAULT;
00493
00494 imgSetStatus(img, STATUS_INVALIDHEADER);
00495
00496
00497 dimNr=h->dime.dim[0];
00498 if(dimNr<2) return STATUS_INVALIDHEADER;
00499 dimx=h->dime.dim[1]; dimy=h->dime.dim[2];
00500 if(dimNr>2) {dimz=h->dime.dim[3]; if(dimNr>3) dimt=h->dime.dim[4];}
00501 pxlNr=dimx*dimy*dimz;
00502 if(pxlNr<1) return STATUS_INVALIDHEADER;
00503 img->dimx=dimx; img->dimy=dimy; img->dimz=dimz; img->dimt=dimt;
00504
00505
00506 img->type=IMG_TYPE_IMAGE;
00507 strncpy(img->studyNr, h->hist.patient_id, MAX_STUDYNR_LEN);
00508 strcpy(img->patientName, h->hist.patient_id);
00509 img->sizex=h->dime.pixdim[1];
00510 img->sizey=h->dime.pixdim[2];
00511 img->sizez=h->dime.pixdim[3];
00512
00513 if(h->dime.funused3>1.E-5) img->isotopeHalflife=h->dime.funused3;
00514 if(h->little) img->_fileFormat=IMG_ANA_L; else img->_fileFormat=IMG_ANA;
00515
00516
00517 if(strstr(h->hist.descrip, "Decay corrected.")!=NULL)
00518 img->decayCorrected=1;
00519 else if(strstr(h->hist.descrip, "No decay correction.")!=NULL)
00520 img->decayCorrected=0;
00521 else
00522 img->decayCorrected=1;
00523
00524 imgSetStatus(img, STATUS_OK);
00525 return STATUS_OK;
00526 }
00527
00528
00529
00542 int imgSetAnalyzeHeader(IMG *img, const char *dbname, ANALYZE_DSR *dsr, float fmin, float fmax) {
00543 struct tm *st;
00544 char *cptr;
00545 float g;
00546
00547 if(IMG_TEST) printf("\nimgSetAnalyzeHeader(*img, *dsr)\n");
00548
00549
00550 if(img==NULL) return STATUS_FAULT;
00551 if(img->status!=IMG_STATUS_INITIALIZED && img->status!=IMG_STATUS_OCCUPIED)
00552 return STATUS_FAULT;
00553 imgSetStatus(img, STATUS_FAULT);
00554 if(dsr==NULL) return STATUS_FAULT;
00555
00556
00557 if(img->_fileFormat==IMG_ANA_L) dsr->little=1; else dsr->little=0;
00558
00559 memset(&dsr->hk, 0, sizeof(ANALYZE_HEADER_KEY));
00560 memset(&dsr->dime, 0, sizeof(ANALYZE_HEADER_IMGDIM));
00561 memset(&dsr->hist, 0, sizeof(ANALYZE_HEADER_HISTORY));
00562 dsr->hk.sizeof_hdr=348;
00563 strcpy(dsr->hk.data_type, "");
00564 cptr=strrchr(dbname, '/'); if(cptr==NULL) cptr=strrchr(dbname, '\\');
00565 if(cptr!=NULL) cptr++; if(cptr==NULL) cptr=(char*)dbname;
00566 strncpy(dsr->hk.db_name, cptr, 17);
00567 dsr->hk.extents=16384;
00568 dsr->hk.regular='r';
00569
00570 dsr->dime.dim[0]=4;
00571 dsr->dime.dim[1]=img->dimx;
00572 dsr->dime.dim[2]=img->dimy;
00573 dsr->dime.dim[3]=img->dimz;
00574 dsr->dime.dim[4]=img->dimt;
00575 dsr->dime.datatype=ANALYZE_DT_SIGNED_SHORT;
00576 dsr->dime.bitpix=16;
00577 dsr->dime.pixdim[0]=0.0;
00578 dsr->dime.pixdim[1]=img->sizex;
00579 dsr->dime.pixdim[2]=img->sizey;
00580 dsr->dime.pixdim[3]=img->sizez;
00581 dsr->dime.pixdim[4]=0.0;
00582 dsr->dime.funused1=0.0;
00583
00584 dsr->dime.funused3=img->isotopeHalflife;
00585
00586 if(img->decayCorrected==1) strcpy(dsr->hist.descrip, "Decay corrected.");
00587 else strcpy(dsr->hist.descrip, "No decay correction.");
00588 strncpy(dsr->hist.scannum, img->studyNr, 10);
00589 st=localtime(&img->scanStart);
00590 if(st!=NULL) {
00591 strftime(dsr->hist.exp_date, 10, "%Y-%m-%d", st);
00592 strftime(dsr->hist.exp_time, 10, "%H:%M:%S", st);
00593 } else {
00594 strcpy(dsr->hist.exp_date, "1900-01-01");
00595 strcpy(dsr->hist.exp_time, "00:00:00");
00596 }
00597
00598 if(fmin<fmax) {
00599 dsr->dime.cal_min=fmin; dsr->dime.cal_max=fmax;
00600 } else {
00601 if(img->status==IMG_STATUS_OCCUPIED &&
00602 imgMinMax(img, &dsr->dime.cal_min, &dsr->dime.cal_max)==0) {}
00603 else return STATUS_FAULT;
00604 }
00605 if(fabs(dsr->dime.cal_min) > fabs(dsr->dime.cal_max)) g = fabs(dsr->dime.cal_min);
00606 else g = fabs(dsr->dime.cal_max);
00607
00608
00609 if(g<1E-20) g=1.0; else g=32767./g; dsr->dime.funused1=1.0/g;
00610
00611 dsr->dime.glmin=temp_roundf(fmin*g); dsr->dime.glmax=temp_roundf(fmax*g);
00612
00613
00614
00615 imgSetStatus(img, STATUS_OK);
00616 return STATUS_OK;
00617 }
00618
00619
00620
00629 int imgReadAnalyzeFirstFrame(const char *fname, IMG *img) {
00630 int ret=0;
00631
00632 if(IMG_TEST) printf("\nimgReadAnalyzeFirstFrame(%s, *img)\n", fname);
00633
00634 if(img==NULL) return STATUS_FAULT;
00635 if(img->status!=IMG_STATUS_INITIALIZED) return STATUS_FAULT;
00636 imgSetStatus(img, STATUS_FAULT);
00637 if(fname==NULL) return STATUS_FAULT;
00638
00639
00640 ret=imgReadAnalyzeHeader(fname, img); if(ret) return(ret);
00641 if(IMG_TEST>3) imgInfo(img);
00642
00643
00644 img->dimt=1;
00645 ret=imgAllocate(img, img->dimz, img->dimy, img->dimx, img->dimt);
00646 if(ret) return STATUS_NOMEMORY;
00647
00648
00649 ret=imgReadAnalyzeFrame(fname, 1, img, 0); if(ret) return(ret);
00650
00651
00652 imgSetStatus(img, STATUS_OK);
00653 return STATUS_OK;
00654 }
00655
00656
00657
00674 int imgReadAnalyzeFrame(const char *fname, int frame_to_read, IMG *img, int frame_index) {
00675 FILE *fp;
00676 int ret, pi, xi, yi;
00677 float *fdata=NULL, *fptr;
00678 ANALYZE_DSR dsr;
00679 char datfile[FILENAME_MAX], hdrfile[FILENAME_MAX], siffile[FILENAME_MAX];
00680 SIF sif;
00681
00682
00683 if(IMG_TEST) printf("\nimgReadAnalyzeFrame(%s, %d, *img, %d)\n",
00684 fname, frame_to_read, frame_index);
00685
00686
00687 if(img==NULL) return STATUS_FAULT;
00688 if(img->status!=IMG_STATUS_OCCUPIED) return STATUS_FAULT;
00689 if(fname==NULL) return STATUS_FAULT;
00690 if(frame_index<0 || frame_index>img->dimt-1) return STATUS_FAULT;
00691 if(frame_to_read<1) return STATUS_FAULT;
00692 imgSetStatus(img, STATUS_FAULT);
00693
00694
00695 ret=anaDatabaseExists(fname, hdrfile, datfile, siffile);
00696 if(ret==0) return STATUS_NOFILE;
00697
00698
00699 ret=anaReadHeader(hdrfile, &dsr);
00700 if(ret!=0) {
00701 if(ret==1) return STATUS_FAULT;
00702 else if(ret==2) return STATUS_NOHEADERFILE;
00703 else return STATUS_UNSUPPORTED;
00704 return(STATUS_FAULT);
00705 }
00706
00707
00708 if(IMG_TEST>2) fprintf(stdout, "reading image data %s\n", datfile);
00709 if((fp=fopen(datfile, "rb")) == NULL) return STATUS_NOIMGDATA;
00710
00711
00712 fdata=malloc(img->dimx*img->dimy*img->dimz*sizeof(float));
00713 if(fdata==NULL) {fclose(fp); return STATUS_NOMEMORY;}
00714
00715
00716 fptr=fdata;
00717 ret=anaReadImagedata(fp, &dsr, frame_to_read, fptr);
00718 fclose(fp);
00719 if(ret==3) {free(fdata); return STATUS_NOMATRIX;}
00720 if(ret!=0) {free(fdata); return STATUS_UNSUPPORTED;}
00721
00722
00723 fptr=fdata;
00724 if(anaFlipping()==0) {
00725 for(pi=0; pi<img->dimz; pi++)
00726 for(yi=img->dimy-1; yi>=0; yi--)
00727 for(xi=img->dimx-1; xi>=0; xi--)
00728 img->m[pi][yi][xi][frame_index]=*fptr++;
00729 } else {
00730 for(pi=img->dimz-1; pi>=0; pi--)
00731 for(yi=img->dimy-1; yi>=0; yi--)
00732 for(xi=img->dimx-1; xi>=0; xi--)
00733 img->m[pi][yi][xi][frame_index]=*fptr++;
00734 }
00735 free(fdata);
00736
00737 imgSetStatus(img, STATUS_OK);
00738
00739
00740
00741
00742 sifInit(&sif);
00743 if(sifRead(siffile, &sif)!=0) return STATUS_OK;
00744
00745 if(sif.frameNr>=frame_to_read) {
00746 img->start[frame_index]=sif.x1[frame_to_read-1];
00747 img->end[frame_index]=sif.x2[frame_to_read-1];
00748 img->mid[frame_index]=0.5*(img->start[frame_index]+img->end[frame_index]);
00749 img->prompts[frame_index]=sif.prompts[frame_to_read-1];
00750 img->randoms[frame_index]=sif.randoms[frame_to_read-1];
00751 }
00752 sifEmpty(&sif);
00753
00754 return STATUS_OK;
00755 }
00756
00757
00758
00781 int imgWriteAnalyzeFrame(const char *dbname, int frame_to_write, IMG *img, int frame_index, float fmin, float fmax) {
00782 IMG test_img;
00783 int ret=0, pxlNr, zi, xi, yi, little;
00784 FILE *fp;
00785 short int *sdata=NULL, *sptr;
00786 char datfile[FILENAME_MAX], hdrfile[FILENAME_MAX], siffile[FILENAME_MAX];
00787 ANALYZE_DSR dsr;
00788 float scale_factor=1.0;
00789
00790
00791 if(IMG_TEST) printf("\nimgWriteAnalyzeFrame(%s, %d, *img, %d, %g, %g)\n",
00792 dbname, frame_to_write, frame_index, fmin, fmax);
00793
00794
00795
00796
00797 if(dbname==NULL) return STATUS_FAULT;
00798 if(img==NULL) return STATUS_FAULT;
00799 if(img->status!=IMG_STATUS_OCCUPIED) return STATUS_FAULT;
00800 if(frame_to_write<0) return STATUS_FAULT;
00801 if(frame_index<0 || frame_index>=img->dimt) return STATUS_FAULT;
00802 if(img->_fileFormat!=IMG_ANA_L && img->_fileFormat!=IMG_ANA) return STATUS_FAULT;
00803
00804
00805
00806
00807
00808
00809
00810 imgInit(&test_img);
00811 if(anaDatabaseExists(dbname, hdrfile, datfile, siffile)==0) {
00812
00813
00814 sprintf(hdrfile, "%s.hdr", dbname);
00815 sprintf(datfile, "%s.img", dbname);
00816 sprintf(siffile, "%s.sif", dbname);
00817
00818
00819 imgSetAnalyzeHeader(img, dbname, &dsr, fmin, fmax);
00820 if(frame_to_write==0) frame_to_write=1;
00821 dsr.dime.dim[4]=frame_to_write;
00822 scale_factor=dsr.dime.funused1;
00823 if(fabs(scale_factor)>1.0E-20) scale_factor=1.0/scale_factor;
00824
00825
00826 ret=anaWriteHeader(hdrfile, &dsr); if(ret && IMG_TEST) printf("anaWriteHeader() := %d\n", ret);
00827 if(ret) return STATUS_CANTWRITEHEADERFILE;
00828
00829
00830 if(access(datfile, 0) != -1) remove(datfile);
00831
00832 } else {
00833
00834
00835 ret=imgReadAnalyzeHeader(dbname, &test_img); if(ret!=0) return ret;
00836
00837 if(img->_fileFormat!=test_img._fileFormat || img->type!=test_img.type)
00838 return STATUS_WRONGFILETYPE;
00839
00840 if(img->dimz!=test_img.dimz || img->dimx!=test_img.dimx ||
00841 img->dimy!=test_img.dimy)
00842 return STATUS_VARMATSIZE;
00843 imgEmpty(&test_img);
00844
00845
00846
00847 if((ret=anaReadHeader(hdrfile, &dsr))!=0) return STATUS_NOMAINHEADER;
00848 scale_factor=1.0/dsr.dime.funused1;
00849 if(frame_to_write==0) frame_to_write=dsr.dime.dim[4]+1;
00850 if(dsr.dime.dim[4]<frame_to_write) {
00851 if(dsr.dime.dim[4]+1<frame_to_write) return STATUS_MISSINGMATRIX;
00852 dsr.dime.dim[4]=frame_to_write;
00853 }
00854 if((ret=anaWriteHeader(hdrfile, &dsr))!=0) return STATUS_NOWRITEPERM;
00855 }
00856 if(IMG_TEST>2) {
00857 printf("frame_to_write := %d\n", frame_to_write);
00858 printf("hdrfile := %s\n", hdrfile);
00859 printf("datfile := %s\n", datfile);
00860 printf("siffile := %s\n", siffile);
00861 }
00862
00863
00864 pxlNr=img->dimx*img->dimy;
00865 sdata=(short int*)malloc(pxlNr*sizeof(short int));
00866 if(sdata==NULL) return STATUS_NOMEMORY;
00867
00868
00869 if(frame_to_write==1) fp=fopen(datfile, "wb"); else fp=fopen(datfile, "r+b");
00870 if(fp==NULL) {free(sdata); return STATUS_CANTWRITEIMGFILE;}
00871
00872 if(fseek(fp, (frame_to_write-1)*pxlNr*img->dimz*sizeof(short int), SEEK_SET)!=0) {
00873 free(sdata); fclose(fp); return STATUS_MISSINGMATRIX;}
00874 little=little_endian();
00875
00876 if(anaFlipping()==0) {
00877 for(zi=0; zi<img->dimz; zi++) {
00878 sptr=sdata;
00879 for(yi=img->dimy-1; yi>=0; yi--) for(xi=img->dimx-1; xi>=0; xi--) {
00880 *sptr=temp_roundf(scale_factor*img->m[zi][yi][xi][frame_index]); sptr++;
00881 }
00882
00883 sptr=sdata; if(little!=dsr.little) swabip(sptr, pxlNr*sizeof(short int));
00884
00885 sptr=sdata;
00886 if(fwrite(sptr, sizeof(short int), pxlNr, fp) != pxlNr) {
00887 free(sdata); fclose(fp); return STATUS_CANTWRITEIMGFILE;
00888 }
00889 }
00890 } else {
00891 for(zi=img->dimz-1; zi>=0; zi--) {
00892 sptr=sdata;
00893 for(yi=img->dimy-1; yi>=0; yi--) for(xi=img->dimx-1; xi>=0; xi--) {
00894 *sptr=temp_roundf(scale_factor*img->m[zi][yi][xi][frame_index]); sptr++;
00895 }
00896
00897 sptr=sdata; if(little!=dsr.little) swabip(sptr, pxlNr*sizeof(short int));
00898
00899 sptr=sdata;
00900 if(fwrite(sptr, sizeof(short int), pxlNr, fp) != pxlNr) {
00901 free(sdata); fclose(fp); return STATUS_CANTWRITEIMGFILE;
00902 }
00903 }
00904 }
00905 free(sdata);
00906 fclose(fp);
00907
00908 return STATUS_OK;
00909 }
00910
00911
00912