/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * filename: m-anlz.c * * * * UTIL C-source: Medical Image Conversion Utility * * * * purpose : read and write ANALYZE files * * * * project : (X)MedCon by Erik Nolf * * * * Functions : MdcCheckANLZ() - Check ANALYZE format * * MdcReadANLZ() - Read ANALYZE file * * MdcWriteANLZ() - Write ANALYZE file * * MdcWriteHeaderKey() - Write Header Key to file * * MdcWriteImageDimension() - Write Image Dimension to file * * MdcWriteDataHistory() - Write Data History to file * * MdcWriteImagesData() - Write the images to file * * MdcGetSpmOpt() - Get specific SPM options * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */ /* $Id: m-anlz.c,v 1.80 2007/05/21 20:16:10 enlf Exp $ */ /* Copyright (C) 1997-2007 by Erik Nolf This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program; if not, write to the Free Software Foundation, Inc., 59 Place - Suite 330, Boston, MA 02111-1307, USA. */ /**************************************************************************** H E A D E R S ****************************************************************************/ #include "m-depend.h" #include #ifdef HAVE_STDLIB_H #include #endif #ifdef HAVE_STRING_H #include #endif #ifdef HAVE_STRINGS_H #ifndef _WIN32 #include #endif #endif #ifdef HAVE_UNISTD_H #include #endif #include "medcon.h" /**************************************************************************** D E F I N E S *****************************************************************************/ static Int8 INIT_SPMOPT = MDC_YES; static MDC_SPMOPT spmopt; #define MDC_ALWAYS_SET_4D 1 /* 0/1 disable/enable always set 4 dims */ /**************************************************************************** F U N C T I O N S *****************************************************************************/ int MdcCheckANLZ(FILEINFO *fi) { MDC_ANLZ_HEADER_KEY hk; int check=2, FORMAT=MDC_FRMT_NONE; if (fread((char *)&hk,1,MDC_ANLZ_HK_SIZE,fi->ifp) != MDC_ANLZ_HK_SIZE) return(MDC_BAD_READ); MDC_FILE_ENDIAN = MDC_HOST_ENDIAN; while (check--) { if ( (hk.sizeof_hdr==348 || hk.sizeof_hdr==148 || hk.sizeof_hdr==228 || hk.sizeof_hdr==384) && (hk.regular == MDC_ANLZ_SIG ) ) { FORMAT = MDC_FRMT_ANLZ; break; } MDC_FILE_ENDIAN = !MDC_HOST_ENDIAN; MdcSWAP(hk.sizeof_hdr); } return(FORMAT); } const char *MdcReadANLZ(FILEINFO *fi) { MDC_SPMOPT *opt = &spmopt; FILE *fp=fi->ifp; MDC_ANLZ_HEADER_KEY hk; MDC_ANLZ_IMAGE_DIMS imd; MDC_ANLZ_DATA_HIST dh; IMG_DATA *id=NULL; DYNAMIC_DATA *dd=NULL; Uint32 bytes, i, plane, f, number; Uint8 *img8=NULL; Int8 WAS_COMPRESSED = MDC_NO; char *origpath=NULL; const char *err=NULL; if (MDC_FILE_STDIN == MDC_YES) return("ANLZ File input from stdin unsupported"); if (MDC_PROGRESS) MdcProgress(MDC_PROGRESS_BEGIN,0.,"Reading Analyze:"); if (MDC_VERBOSE) MdcPrntMesg("ANLZ Reading <%s> ...",fi->ifname); /* get endian of the file in MDC_FILE_ENDIAN */ i=MdcCheckANLZ(fi); fseek(fp,0,SEEK_SET); if (i != MDC_FRMT_ANLZ) { if (MDC_FALLBACK_FRMT == MDC_FRMT_ANLZ) { /* set host endian to try analyze */ MDC_FILE_ENDIAN = MDC_HOST_ENDIAN; }else{ /* bail out */ return("ANLZ Endian check failed"); } } memset(&hk,0,MDC_ANLZ_HK_SIZE); memset(&imd,0,MDC_ANLZ_IMD_SIZE); memset(&dh,0,MDC_ANLZ_DH_SIZE); /* put some default we use */ fi->reconstructed = MDC_YES; fi->acquisition_type = MDC_ACQUISITION_TOMO; dh.orient=(char)0xff; if (fread((char *)&hk,1,MDC_ANLZ_HK_SIZE,fp) != MDC_ANLZ_HK_SIZE) return("ANLZ Bad read HeadKey struct"); fi->endian=MDC_FILE_ENDIAN; MdcSWAP(hk.sizeof_hdr); MdcSWAP(hk.extents); MdcSWAP(hk.session_error); if (MDC_INFO) { MdcPrntScrn("\nMDC_ANLZ_HEADER_KEY (%d bytes)\n",MDC_ANLZ_HK_SIZE); MdcPrintLine('-',MDC_HALF_LENGTH); MdcPrntScrn("sizeof_hdr : %d\n",hk.sizeof_hdr); strncpy(mdcbufr,hk.data_type,10); mdcbufr[10]='\0'; MdcPrntScrn("data_type : "); MdcPrintStr(mdcbufr); strncpy(mdcbufr,hk.db_name,18); mdcbufr[18]='\0'; MdcPrntScrn("db_name : "); MdcPrintStr(mdcbufr); MdcPrntScrn("extents : %d\n",hk.extents); MdcPrntScrn("session_error : %hd\n",hk.session_error); MdcPrntScrn("regular : "); MdcPrintChar(hk.regular); MdcPrntScrn("\n"); MdcPrntScrn("hkey_un0 : "); MdcPrintChar(hk.hkey_un0); MdcPrntScrn("\n"); } if (MDC_INFO) { MdcPrntScrn("\nIMAGE_DIMENSION (%d bytes)\n",MDC_ANLZ_IMD_SIZE); MdcPrintLine('-',MDC_HALF_LENGTH); } if (fread((char *)&imd,1,MDC_ANLZ_IMD_SIZE,fp)!=MDC_ANLZ_IMD_SIZE) return("ANLZ Bad read ImageDimensions struct"); for (i=0; i (MDC_ANLZ_MAX_DIMS - 1))) { if (MDC_FALLBACK_FRMT == MDC_FRMT_ANLZ) { /* force reading, set to an acceptable value */ /* hope unused dims were initialized to zero */ for (i = 3; i < MDC_ANLZ_MAX_DIMS; i++) if (imd.dim[i] <= 0) break; imd.dim[0] = i-1; MdcPrntWarn("ANLZ Bad header value in dim[0] dimension"); }else{ /* bail out safely */ return("ANLZ Bad header value in dim[0] dimension"); } } if (MDC_INFO) { for (i=0; i - ANLZ Truncated header",fi->ifname); } memcpy(&opt->origin_x,&dh.originator[0],2); MdcSWAP(opt->origin_x); memcpy(&opt->origin_y,&dh.originator[2],2); MdcSWAP(opt->origin_y); memcpy(&opt->origin_z,&dh.originator[4],2); MdcSWAP(opt->origin_z); MdcSWAP(dh.views); MdcSWAP(dh.vols_added); MdcSWAP(dh.start_field); MdcSWAP(dh.field_skip); MdcSWAP(dh.omax); MdcSWAP(dh.omin); MdcSWAP(dh.smax); MdcSWAP(dh.smin); if (MDC_INFO) { strncpy(mdcbufr,dh.descrip,80); mdcbufr[80]='\0'; MdcPrntScrn("description : "); MdcPrintStr(mdcbufr); strncpy(mdcbufr,dh.aux_file,24); mdcbufr[24]='\0'; MdcPrntScrn("aux_file : "); MdcPrintStr(mdcbufr); MdcPrntScrn("orient : "); switch (dh.orient) { case MDC_ANLZ_TRANS_UNFLIPPED: MdcPrntScrn("transverse unflipped"); break; case MDC_ANLZ_CORON_UNFLIPPED: MdcPrntScrn("coronal unflipped"); break; case MDC_ANLZ_SAGIT_UNFLIPPED: MdcPrntScrn("sagittal unflipped"); break; case MDC_ANLZ_TRANS_FLIPPED : MdcPrntScrn("transverse flipped"); break; case MDC_ANLZ_CORON_FLIPPED : MdcPrntScrn("coronal flipped"); break; case MDC_ANLZ_SAGIT_FLIPPED : MdcPrntScrn("sagittal flipped"); break; default: MdcPrntScrn("Unknown"); } MdcPrntScrn("\n"); strncpy(mdcbufr,dh.originator,10); mdcbufr[10]='\0'; MdcPrntScrn("originator : "); MdcPrintStr(mdcbufr); strncpy(mdcbufr,dh.generated,10); mdcbufr[10]='\0'; MdcPrntScrn("generated : "); MdcPrintStr(mdcbufr); strncpy(mdcbufr,dh.scannum,10); mdcbufr[10]='\0'; MdcPrntScrn("scannum : "); MdcPrintStr(mdcbufr); strncpy(mdcbufr,dh.patient_id,10); mdcbufr[10]='\0'; MdcPrntScrn("patient_id : "); MdcPrintStr(mdcbufr); strncpy(mdcbufr,dh.exp_date,10); mdcbufr[10]='\0'; MdcPrntScrn("exp_date : "); MdcPrintStr(mdcbufr); strncpy(mdcbufr,dh.exp_time,10); mdcbufr[10]='\0'; MdcPrntScrn("exp_time : "); MdcPrintStr(mdcbufr); MdcPrntScrn("views : %d\n",dh.views); MdcPrntScrn("vols_added : %d\n",dh.vols_added); MdcPrntScrn("start_field : %d\n",dh.start_field); MdcPrntScrn("omax : %d\n",dh.omax); MdcPrntScrn("omin : %d\n",dh.omin); MdcPrntScrn("smax : %d\n",dh.smax); MdcPrntScrn("smin : %d\n",dh.smin); } if (MDC_INFO) { MdcPrntScrn("\n"); MdcPrintLine('=',MDC_FULL_LENGTH); MdcPrntScrn("SPM - HEADER INTERPRETATION\n"); MdcPrintLine('-',MDC_HALF_LENGTH); MdcPrntScrn("image {x} : %hd\n",imd.dim[1]); MdcPrntScrn("image {y} : %hd\n",imd.dim[2]); MdcPrntScrn("image {z} : %hd\n",imd.dim[3]); MdcPrntScrn("voxel {x} : %+e\n",imd.pixdim[1]); MdcPrntScrn("voxel {y} : %+e\n",imd.pixdim[2]); MdcPrntScrn("voxel {z} : %+e\n",imd.pixdim[3]); MdcPrntScrn("scaling : %+e\n",imd.spm_pix_rescale); MdcPrntScrn("data type : %hd\n",imd.datatype); MdcPrntScrn("offset : %+e\n",imd.avw_vox_offset); MdcPrntScrn("origin : %hd %hd %hd\n",opt->origin_x ,opt->origin_y ,opt->origin_z); MdcPrntScrn("description : "); MdcPrintStr(dh.descrip); MdcPrintLine('=',MDC_FULL_LENGTH); } /* save the offset, valid for AVW / SPM / MRIcro Analyze files */ opt->offset = imd.avw_vox_offset; /* update our FILEINFO structure */ MdcStringCopy(fi->study_descr,dh.descrip,80); MdcStringCopy(fi->patient_id,dh.patient_id,10); MdcStringCopy(fi->study_id,dh.scannum,10); if (MDC_ECHO_ALIAS == MDC_YES) { MdcEchoAliasName(fi); return(NULL); } memcpy(fi->dim,imd.dim,sizeof(imd.dim)); memcpy(fi->pixdim,imd.pixdim,sizeof(imd.pixdim)); fi->mwidth = (Uint32) imd.dim[1]; fi->mheight = (Uint32) imd.dim[2]; for ( number=1, i=3; i<=imd.dim[0]; i++) number*=imd.dim[i]; if (number == 0) return("ANLZ No valid images specified"); fi->bits = imd.bitpix; switch (imd.datatype) { case MDC_ANLZ_DT_BINARY : fi->type=BIT1; fi->bits=8; break; case MDC_ANLZ_DT_UNSIGNED_CHAR: fi->type=BIT8_U; fi->bits=8; break; case MDC_ANLZ_DT_SIGNED_SHORT : fi->type=BIT16_S; fi->bits=16; break; case MDC_ANLZ_DT_SIGNED_INT : fi->type=BIT32_S; fi->bits=32; break; case MDC_ANLZ_DT_FLOAT : fi->type=FLT32; fi->bits=32; break; case MDC_ANLZ_DT_COMPLEX : return("ANLZ Datatype `complex' unsupported"); break; case MDC_ANLZ_DT_DOUBLE : fi->type=FLT64; fi->bits=64; break; case MDC_ANLZ_DT_RGB : fi->type=COLRGB; fi->bits=24; fi->map=MDC_MAP_PRESENT; break; case MDC_ANLZ_DT_ALL : return("ANLZ Datatype `All' unsupported"); break; default : switch (fi->bits) { case 1: fi->type=BIT1; break; case 8: fi->type=BIT8_U; break; case 16: fi->type=BIT16_S; break; case 32: fi->type=BIT32_S; break; /* could be FLT32 as well */ default: MdcPrntWarn("ANLZ Unknown datatype"); } } /* preserve original path */ MdcMergePath(fi->ipath,fi->idir,fi->ifname); if ((origpath=malloc(strlen(fi->ipath) + 1)) == NULL) return("ANLZ Couldn't allocate original path"); strcpy(origpath,fi->ipath); MdcSplitPath(fi->ipath,fi->idir,fi->ifname); /* read the image file */ MdcCloseFile(fi->ifp); MdcMergePath(fi->ipath,fi->idir,fi->ifname); MdcSetExt(fi->ipath,"img"); /* check for compressed image file */ if (MdcFileExists(fi->ipath) == MDC_NO) { MdcAddCompressionExt(fi->compression,fi->ipath); if (MdcDecompressFile(fi->ipath) != MDC_OK) { MdcFree(origpath); return("ANLZ Decompression image file failed"); } WAS_COMPRESSED = MDC_YES; } if ( (fi->ifp=fopen(fi->ipath,"rb")) == NULL ) { MdcFree(origpath); return("ANLZ Couldn't open image file"); } if (WAS_COMPRESSED == MDC_YES) { unlink(fi->ipath); if (MDC_PROGRESS) MdcProgress(MDC_PROGRESS_BEGIN,0.,"Reading Analyze:"); } MdcSplitPath(fi->ipath,fi->idir,fi->ifname); if (MDC_ANLZ_SPM == MDC_YES) { /* interpret offset value from the header but we keep */ /* our precautions ... badly initialized headers */ long offsetu = (long)opt->offset; if ((float)offsetu == opt->offset) fseek(fi->ifp,offsetu,SEEK_SET); } if (!MdcGetStructID(fi,number)) { MdcFree(origpath); return("ANLZ Bad malloc IMG_DATA structs"); } /* attempt to fill in orientation information */ switch (dh.orient) { /* flipped, what's the meaning of flipped here ? */ case MDC_ANLZ_TRANS_UNFLIPPED: fi->pat_slice_orient=MDC_SUPINE_HEADFIRST_TRANSAXIAL; break; case MDC_ANLZ_CORON_UNFLIPPED: fi->pat_slice_orient=MDC_SUPINE_HEADFIRST_CORONAL; break; case MDC_ANLZ_SAGIT_UNFLIPPED: fi->pat_slice_orient=MDC_SUPINE_HEADFIRST_SAGITTAL; break; case MDC_ANLZ_TRANS_FLIPPED: fi->pat_slice_orient=MDC_SUPINE_HEADFIRST_TRANSAXIAL; break; case MDC_ANLZ_CORON_FLIPPED: fi->pat_slice_orient=MDC_SUPINE_HEADFIRST_CORONAL; break; case MDC_ANLZ_SAGIT_FLIPPED: fi->pat_slice_orient=MDC_SUPINE_HEADFIRST_SAGITTAL; break; } strcpy(fi->pat_pos,MdcGetStrPatPos(fi->pat_slice_orient)); strcpy(fi->pat_orient,MdcGetStrPatOrient(fi->pat_slice_orient)); for ( i=0; i < fi->number; i++) { if (MDC_PROGRESS) MdcProgress(MDC_PROGRESS_INCR,1./(float)fi->number,NULL); plane = i % fi->dim[3]; id = &fi->image[i]; id->width = fi->mwidth; id->height = fi->mheight; id->bits = fi->bits; id->type = fi->type; if (MDC_ANLZ_SPM) { /* consider the scaling factor */ if (imd.spm_pix_rescale > 0.0 ) id->quant_scale = imd.spm_pix_rescale; } if (fi->pixdim[0] == 3.0 ) { id->pixel_xsize = fi->pixdim[1]; id->pixel_ysize = fi->pixdim[2]; id->slice_width = fi->pixdim[3]; }else if (fi->pixdim[0] == 4.0 ) { id->pixel_xsize = fi->pixdim[1]; id->pixel_ysize = fi->pixdim[2]; id->slice_width = fi->pixdim[3]; }else if ( (fi->pixdim[1] > 0.0) && (fi->pixdim[2] > 0.0) && (fi->pixdim[3] > 0.0) ) { /* we will try it anyway */ /* some don't fill in pixdim[0] */ /* for example PMOD (11-Apr-2000) */ id->pixel_xsize = fi->pixdim[1]; id->pixel_ysize = fi->pixdim[2]; id->slice_width = fi->pixdim[3]; fi->pixdim[0] = 3.0; }else { id->pixel_xsize = 1.0; id->pixel_ysize = 1.0; id->slice_width = 1.0; } id->slice_spacing = id->slice_width; MdcFillImgPos(fi,i,plane,0.0); MdcFillImgOrient(fi,i); bytes = MdcPixels2Bytes(fi->mwidth*fi->mheight*fi->bits); if ( (id->buf=MdcGetImgBuffer(bytes)) == NULL ) { MdcFree(img8); MdcFree(origpath); return("ANLZ Bad malloc image buffer"); } if (img8 != NULL) { /* image from buffer */ memcpy(id->buf,img8+i*bytes,bytes); }else{ /* image from file */ if (fread(id->buf,1,bytes,fi->ifp) != bytes ) { err=MdcHandleTruncated(fi, i+1,MDC_YES); if (err != NULL) { MdcFree(origpath); return(err); } } if (fi->truncated) break; } } MdcFree(img8); MdcCloseFile(fi->ifp); /* check some final FILEINFO entries */ if (fi->dim[4] > 1) { fi->acquisition_type = MDC_ACQUISITION_DYNAMIC; /* fill in dynamic data struct */ if (!MdcGetStructDD(fi,(unsigned)fi->dim[4])) return("ANLZ Couldn't malloc DYNAMIC_DATA structs"); for (f=0; f < fi->dynnr; f++) { dd = &fi->dyndata[f]; dd->nr_of_slices = fi->dim[3]; dd->time_frame_duration = fi->pixdim[4]; dd->time_frame_start = f * dd->time_frame_duration; } } /* restore original filename */ strcpy(fi->ipath,origpath); MdcSplitPath(fi->ipath,fi->idir,fi->ifname); MdcFree(origpath); if (fi->truncated) return("ANLZ Truncated image file"); return NULL; } int MdcWriteHeaderKey(FILEINFO *fi) { MDC_ANLZ_HEADER_KEY hk; char *p = NULL; memset(&hk,0,MDC_ANLZ_HK_SIZE); hk.sizeof_hdr = MDC_ANLZ_HK_SIZE + MDC_ANLZ_IMD_SIZE + MDC_ANLZ_DH_SIZE; sprintf(hk.data_type,"dsr"); MdcSplitPath(fi->opath,fi->odir,fi->ofname); p = strrchr(fi->ofname,'.'); if (p != NULL) *p = '\0'; /* remove extension */ sprintf(hk.db_name,"%.18s",fi->ofname); if (p != NULL) *p = '.'; /* add extension */ MdcMergePath(fi->opath,fi->odir,fi->ofname); hk.extents=16384; hk.session_error=0; hk.regular='r'; MdcSWAP(hk.sizeof_hdr); MdcSWAP(hk.extents); MdcSWAP(hk.session_error); fwrite((char *)&hk,1,MDC_ANLZ_HK_SIZE,fi->ofp); if (ferror(fi->ofp)) return(MDC_NO); return(MDC_YES); } int MdcWriteImageDimension(FILEINFO *fi, MDC_SPMOPT *opt) { MDC_ANLZ_IMAGE_DIMS imd; float glmax, glmin; int i; memset(&imd,0,MDC_ANLZ_IMD_SIZE); strcpy(imd.avw_vox_units,"mm"); for (i=0; i <= fi->dim[0]; i++) imd.dim[i] = fi->dim[i]; for (i=0; i <= fi->pixdim[0]; i++) imd.pixdim[i] = fi->pixdim[i]; #if MDC_ALWAYS_SET_4D /* set dummy 4th dimension (time) */ if (imd.dim[0] == 3) { imd.dim[0] = 4; imd.dim[4] = 1; } if (imd.pixdim[0] == 3.) { imd.pixdim[0] = 4.; imd.pixdim[4] = 0.; } #endif #ifdef MDC_USE_SLICE_SPACING if (fi->number > 1) imd.pixdim[3] = fi->image[0].slice_spacing; #endif imd.dim[1] = (Int16) fi->mwidth; imd.dim[2] = (Int16) fi->mheight; if (fi->map == MDC_MAP_PRESENT) { /* colored */ imd.datatype = MDC_ANLZ_DT_RGB; imd.bitpix = 24; }else{ /* grayscale */ if (MDC_FORCE_INT != MDC_NO) { switch (MDC_FORCE_INT) { case BIT8_U : imd.datatype = MDC_ANLZ_DT_UNSIGNED_CHAR; imd.bitpix = 8; break; case BIT16_S: imd.datatype = MDC_ANLZ_DT_SIGNED_SHORT; imd.bitpix = 16; break; default : imd.datatype = MDC_ANLZ_DT_SIGNED_SHORT; imd.bitpix = 16; } }else if (!(MDC_QUANTIFY || MDC_CALIBRATE)) { if ( fi->diff_type ) { imd.datatype = MDC_ANLZ_DT_SIGNED_SHORT; imd.bitpix = 16; }else{ switch ( fi->type ) { case BIT8_U: case BIT8_S: imd.datatype = MDC_ANLZ_DT_UNSIGNED_CHAR; imd.bitpix = 8; break; case BIT16_U: case BIT16_S: imd.datatype = MDC_ANLZ_DT_SIGNED_SHORT; imd.bitpix = 16; break; #ifdef HAVE_8BYTE_INT case BIT64_U: case BIT64_S: #endif case BIT32_U: case BIT32_S: imd.datatype = MDC_ANLZ_DT_SIGNED_INT; imd.bitpix = 32; break; case FLT32: imd.datatype = MDC_ANLZ_DT_FLOAT; imd.bitpix = 32; break; case FLT64: imd.datatype = MDC_ANLZ_DT_DOUBLE; imd.bitpix = 64; break; } } }else{ if (MDC_ANLZ_SPM == MDC_YES) { /* BIT16_S with scaling factor */ imd.datatype = MDC_ANLZ_DT_SIGNED_SHORT; imd.bitpix = 16; }else{ imd.datatype = MDC_ANLZ_DT_FLOAT; imd.bitpix = 32; } } } /* find and set max/min values */ for (i = 0; i < fi->number; i++) { IMG_DATA *id = &fi->image[i]; if (id->rescaled == MDC_YES) { if (i == 0) { /* init values */ glmax = id->rescaled_max; glmin = id->rescaled_min; }else{ /* get max/min */ glmax = (id->rescaled_max > glmax) ? id->rescaled_max : glmax; glmin = (id->rescaled_min < glmin) ? id->rescaled_min : glmin; } }else{ if (i == 0) { /* init values */ glmax = id->max; glmin = id->min; }else{ /* get max/min */ glmax = (id->max > glmax) ? id->max : glmax; glmin = (id->min < glmin) ? id->min : glmin; } } } imd.glmax = (Int32) glmax; imd.glmin = (Int32) glmin; imd.avw_cal_max = fi->qglmax; imd.avw_cal_min = fi->qglmin; /* thinking about SPM */ if (imd.pixdim[0] <= 0.0 || imd.pixdim[0] >= (float)MDC_ANLZ_MAX_DIMS) { imd.pixdim[0]=3.; imd.pixdim[1]=1.; imd.pixdim[2]=1.; imd.pixdim[3]=1.; } if (opt != NULL) imd.avw_vox_offset = opt->offset; if (MDC_ANLZ_SPM == MDC_YES) { /* the scaling factor */ if (fi->image[0].rescaled) imd.spm_pix_rescale=(float)fi->image[0].rescaled_fctr; /* did rescale over all images -> all images same factor */ }else{ imd.spm_pix_rescale=1.; } /* swap the data if necessary */ for (i=0; iofp); if (ferror(fi->ofp)) return(MDC_NO); return(MDC_YES); } int MdcWriteDataHistory(FILEINFO *fi, MDC_SPMOPT *opt) { MDC_ANLZ_DATA_HIST dh; memset(&dh,0,MDC_ANLZ_DH_SIZE); sprintf(dh.descrip,"%.35s",fi->study_descr); sprintf(dh.scannum,"%.10s",fi->study_id); sprintf(dh.patient_id,"%.10s",fi->patient_id); sprintf(dh.generated,"%.10s",MDC_PRGR); switch (fi->pat_slice_orient) { case MDC_SUPINE_HEADFIRST_TRANSAXIAL : case MDC_PRONE_HEADFIRST_TRANSAXIAL : case MDC_DECUBITUS_RIGHT_HEADFIRST_TRANSAXIAL: case MDC_DECUBITUS_LEFT_HEADFIRST_TRANSAXIAL : case MDC_SUPINE_FEETFIRST_TRANSAXIAL : case MDC_PRONE_FEETFIRST_TRANSAXIAL : case MDC_DECUBITUS_RIGHT_FEETFIRST_TRANSAXIAL: case MDC_DECUBITUS_LEFT_FEETFIRST_TRANSAXIAL : dh.orient = MDC_ANLZ_TRANS_UNFLIPPED; break; case MDC_SUPINE_HEADFIRST_CORONAL : case MDC_PRONE_HEADFIRST_CORONAL : case MDC_DECUBITUS_RIGHT_HEADFIRST_CORONAL : case MDC_DECUBITUS_LEFT_HEADFIRST_CORONAL : case MDC_SUPINE_FEETFIRST_CORONAL : case MDC_PRONE_FEETFIRST_CORONAL : case MDC_DECUBITUS_RIGHT_FEETFIRST_CORONAL : case MDC_DECUBITUS_LEFT_FEETFIRST_CORONAL : dh.orient = MDC_ANLZ_CORON_UNFLIPPED; break; case MDC_SUPINE_HEADFIRST_SAGITTAL : case MDC_PRONE_HEADFIRST_SAGITTAL : case MDC_DECUBITUS_RIGHT_HEADFIRST_SAGITTAL : case MDC_DECUBITUS_LEFT_HEADFIRST_SAGITTAL : case MDC_SUPINE_FEETFIRST_SAGITTAL : case MDC_PRONE_FEETFIRST_SAGITTAL : case MDC_DECUBITUS_RIGHT_FEETFIRST_SAGITTAL : case MDC_DECUBITUS_LEFT_FEETFIRST_SAGITTAL : dh.orient = MDC_ANLZ_SAGIT_UNFLIPPED; break; } if (opt != NULL) { MdcSWAP(opt->origin_x); memcpy(&dh.originator[0],&opt->origin_x,2); MdcSWAP(opt->origin_y); memcpy(&dh.originator[2],&opt->origin_y,2); MdcSWAP(opt->origin_z); memcpy(&dh.originator[4],&opt->origin_z,2); } fwrite((char *)&dh,1,MDC_ANLZ_DH_SIZE,fi->ofp); if (ferror(fi->ofp)) return(MDC_NO); return(MDC_YES); } char *MdcWriteImagesData(FILEINFO *fi) { double pval; Uint8 grval; Uint32 i, FREE; Uint32 size, n, nr; Uint16 type; Uint8 *buf, *maxbuf; Int8 saved_norm_over_frames=MDC_NORM_OVER_FRAMES; IMG_DATA *id; for (i=fi->number; i>0; i-- ) { if (MDC_PROGRESS) MdcProgress(MDC_PROGRESS_INCR,1./(float)fi->number,NULL); nr = fi->number - i; /* normal planes */ id = &fi->image[nr]; buf = id->buf; FREE = MDC_NO; type = id->type; if (fi->map != MDC_MAP_PRESENT) { /* grayscale */ if (MDC_FORCE_INT != MDC_NO) { if (MDC_ANLZ_SPM) MDC_NORM_OVER_FRAMES = MDC_NO; switch (MDC_FORCE_INT) { case BIT8_U : buf = MdcGetImgBIT8_U(fi,nr); type=BIT8_U; FREE=MDC_YES; break; case BIT16_S: buf = MdcGetImgBIT16_S(fi,nr); type=BIT16_S; FREE=MDC_YES; break; default : buf = MdcGetImgBIT16_S(fi,nr); type=BIT16_S; FREE=MDC_YES; } if (MDC_ANLZ_SPM) MDC_NORM_OVER_FRAMES = saved_norm_over_frames; }else if (!(MDC_QUANTIFY || MDC_CALIBRATE)) { if ( fi->diff_type ) { switch (id->type) { case BIT16_S: buf = id->buf; type = BIT16_S; FREE=MDC_NO; break; default : buf = MdcGetImgBIT16_S(fi,nr); type = BIT16_S; FREE=MDC_YES; break; } }else{ switch (id->type) { case BIT8_S: buf = MdcGetImgBIT8_U(fi,nr); type=BIT8_U ; FREE=MDC_YES; break; case BIT16_U: buf = MdcGetImgBIT16_S(fi,nr); type=BIT16_S; FREE=MDC_YES; break; case BIT32_U: buf = MdcGetImgBIT32_S(fi,nr); type=BIT32_S; FREE=MDC_YES; break; case BIT64_S: case BIT64_U: buf = MdcGetImgBIT32_S(fi,nr); type=BIT32_S; FREE=MDC_YES; break; } } }else{ if (MDC_ANLZ_SPM == MDC_YES) { /* using the global scale factor <=> normalize over ALL images! */ /* so all images have the same scale factor */ MDC_NORM_OVER_FRAMES=MDC_NO; buf = MdcGetImgBIT16_S(fi,nr); type = BIT16_S; FREE=MDC_YES; MDC_NORM_OVER_FRAMES=saved_norm_over_frames; }else{ buf = MdcGetImgFLT32(fi,nr); type=FLT32; FREE=MDC_YES; } } } if (buf == NULL) return("ANLZ Bad malloc image buffer"); if (fi->diff_size) { maxbuf = MdcGetResizedImage(fi, buf, type, nr); if (maxbuf == NULL) return("ANLZ Bad malloc maxbuf"); if (FREE) MdcFree(buf); FREE = MDC_YES; }else{ maxbuf = buf; } size = fi->mwidth * fi->mheight * MdcType2Bytes(type); if (fi->type == COLRGB) { /* true color */ if (fwrite(maxbuf,1,size,fi->ofp) != size) return("ANLZ Bad write RGB buffer"); }else{ for (n=0; n < size; n += MdcType2Bytes(type)) { /* indexed */ pval = MdcGetDoublePixel((Uint8 *)&maxbuf[n],type); if (fi->map == MDC_MAP_PRESENT) { /* colored */ grval = (Uint8)pval; fwrite(&fi->palette[grval * 3 + 0], 1, 1, fi->ofp); /* red */ fwrite(&fi->palette[grval * 3 + 1], 1, 1, fi->ofp); /* green */ fwrite(&fi->palette[grval * 3 + 2], 1, 1, fi->ofp); /* blue */ if (ferror(fi->ofp)) return("ANLZ Bad write colored pixel"); }else{ /* grayscale */ if (!MdcWriteDoublePixel(pval,type,fi->ofp)) return("ANLZ Bad write image pixel"); } } } if (FREE) MdcFree(maxbuf); if (ferror(fi->ofp)) return("ANLZ Bad writing of images"); } return NULL; } void MdcGetSpmOpt(FILEINFO *fi, MDC_SPMOPT *opt) { if (INIT_SPMOPT == MDC_YES) { opt->origin_x = 0; opt->origin_y = 0; opt->origin_z = 0; opt->offset = 0.; INIT_SPMOPT = MDC_NO; } if (MDC_FILE_STDIN == MDC_YES) return; /* stdin already in use */ MdcPrintLine('-',MDC_FULL_LENGTH); MdcPrntScrn("\tSPM OPTIONS\t\tORIG FILE: %s\n",fi->ifname); MdcPrintLine('-',MDC_FULL_LENGTH); MdcPrntScrn("\n\tThe origin values must be an Int16 value"); MdcPrntScrn("\n\tThere is NO check performed on the input!\n"); MdcPrntScrn("\n\tOrigin X [%d]? ",opt->origin_x); if (!MdcPutDefault(mdcbufr)) opt->origin_x = (Int16)atoi(mdcbufr); MdcPrntScrn("\n\tOrigin Y [%d]? ",opt->origin_y); if (!MdcPutDefault(mdcbufr)) opt->origin_y = (Int16)atoi(mdcbufr); MdcPrntScrn("\n\tOrigin Z [%d]? ",opt->origin_z); if (!MdcPutDefault(mdcbufr)) opt->origin_z = (Int16)atoi(mdcbufr); /* MARK: skip asking about offset */ /* MdcPrntScrn("\n\tOffset [%+e]? ",opt->offset); */ /* if (!MdcPutDefault(mdcbufr)) opt->offset = (float)atof(mdcbufr); */ MdcPrntScrn("\n"); MdcPrintLine('-',MDC_FULL_LENGTH); } const char *MdcWriteANLZ(FILEINFO *fi) { MDC_SPMOPT *opt = &spmopt; char tmpfname[MDC_MAX_PATH + 1]; const char *msg; MDC_FILE_ENDIAN = MDC_WRITE_ENDIAN; /* user wanted to supply some parameters */ if ((MDC_ANLZ_OPTIONS == MDC_YES) && (XMDC_GUI == MDC_NO)) { MdcGetSpmOpt(fi,opt); }else { /* set default origin to image centre of middle slice */ opt->origin_x = (Int16)((fi->dim[1] + 1)/2); opt->origin_y = (Int16)((fi->dim[2] + 1)/2); opt->origin_z = (Int16)((fi->dim[3] + 1)/2); opt->offset = 0.; } /* header and image separate, rescaled stuff very important */ /* so we will write the images first ! */ /* get filename: no longer with truncation */ /* SPM, PMOD etc don't rely on db_name[18]) */ if (XMDC_GUI == MDC_YES) { strcpy(tmpfname,fi->opath); }else{ if (MDC_ALIAS_NAME == MDC_YES) { MdcAliasName(fi,tmpfname); }else{ strcpy(tmpfname,fi->ifname); } MdcDefaultName(fi,MDC_FRMT_ANLZ,fi->ofname,tmpfname); } if (MDC_PROGRESS) MdcProgress(MDC_PROGRESS_BEGIN,0.,"Writing Analyze:"); if (MDC_VERBOSE) MdcPrntMesg("ANLZ Writing <%s> & <.img> ...",fi->ofname); /* writing images */ if (XMDC_GUI == MDC_YES) { fi->ofname[0]='\0'; MdcNewExt(fi->ofname,tmpfname,"img"); }else{ MdcNewName(fi->ofname,tmpfname,"img"); } if (MDC_FILE_STDOUT == MDC_YES) { /* send image data to stdout (1>stdout) */ fi->ofp = stdout; }else{ if (MdcKeepFile(fi->ofname)) return("ANLZ Image file exists!!"); if ( (fi->ofp=fopen(fi->ofname,"wb")) == NULL ) return ("ANLZ Couldn't open image file"); } msg = MdcWriteImagesData(fi); if (msg != NULL) return(msg); MdcCloseFile(fi->ofp); /* writing header with rescaled stuff */ if (XMDC_GUI == MDC_YES) { strcpy(fi->ofname,tmpfname); }else{ MdcDefaultName(fi,MDC_FRMT_ANLZ,fi->ofname,tmpfname); } if (MDC_FILE_STDOUT == MDC_YES) { /* send header to stderr (2>stderr) */ fi->ofp = stderr; }else{ if (MdcKeepFile(fi->ofname)) return("ANLZ Header file exists!!"); if ( (fi->ofp=fopen(fi->ofname,"wb")) == NULL ) return("ANLZ Couldn't open header file"); } if ( !MdcWriteHeaderKey(fi) ) return("ANLZ Bad write HeaderKey struct"); if ( !MdcWriteImageDimension(fi, opt) ) return("ANLZ Bad write ImageDimension struct"); if ( !MdcWriteDataHistory(fi, opt) ) return("ANLZ Bad write DataHistory struct"); MdcCheckQuantitation(fi); MdcCloseFile(fi->ofp); return(NULL); }