/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * filename: m-rslice.c * * * * UTIL C-source: Medical Image Conversion Utility * * * * purpose : reslice in different projections * * * * project : (X)MedCon by Erik Nolf * * * * Functions : MdcGetImageProjection() - Retrieve image projection * * MdcGetNewPatSliceOrient() - Get new patient slice orient * * MdcCheckReslice() - Check before reslicing * * MdcResliceImages() - Reslice image (tra,cor,sag) * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */ /* $Id: m-rslice.c,v 1.28 2007/05/21 20:16:14 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 #include "medcon.h" /**************************************************************************** D E F I N E S ****************************************************************************/ #define MDC_VOLUME_REQUIRED 10 /* minimum slices considered viewable */ /**************************************************************************** F U N C T I O N S ****************************************************************************/ Int8 MdcGetSliceProjection(FILEINFO *cur) { Int8 slice_projection = cur->slice_projection; if (cur->slice_projection == MDC_UNKNOWN) { switch (cur->pat_slice_orient) { case MDC_SUPINE_HEADFIRST_TRANSAXIAL : case MDC_SUPINE_FEETFIRST_TRANSAXIAL : case MDC_PRONE_HEADFIRST_TRANSAXIAL : case MDC_PRONE_FEETFIRST_TRANSAXIAL : case MDC_DECUBITUS_RIGHT_HEADFIRST_TRANSAXIAL : case MDC_DECUBITUS_RIGHT_FEETFIRST_TRANSAXIAL : case MDC_DECUBITUS_LEFT_HEADFIRST_TRANSAXIAL : case MDC_DECUBITUS_LEFT_FEETFIRST_TRANSAXIAL : slice_projection = MDC_TRANSAXIAL; break; case MDC_SUPINE_HEADFIRST_SAGITTAL : case MDC_SUPINE_FEETFIRST_SAGITTAL : case MDC_PRONE_HEADFIRST_SAGITTAL : case MDC_PRONE_FEETFIRST_SAGITTAL : case MDC_DECUBITUS_RIGHT_HEADFIRST_SAGITTAL : case MDC_DECUBITUS_RIGHT_FEETFIRST_SAGITTAL : case MDC_DECUBITUS_LEFT_HEADFIRST_SAGITTAL : case MDC_DECUBITUS_LEFT_FEETFIRST_SAGITTAL : slice_projection = MDC_SAGITTAL; break; case MDC_SUPINE_HEADFIRST_CORONAL : case MDC_SUPINE_FEETFIRST_CORONAL : case MDC_PRONE_HEADFIRST_CORONAL : case MDC_PRONE_FEETFIRST_CORONAL : case MDC_DECUBITUS_RIGHT_HEADFIRST_CORONAL : case MDC_DECUBITUS_RIGHT_FEETFIRST_CORONAL : case MDC_DECUBITUS_LEFT_HEADFIRST_CORONAL : case MDC_DECUBITUS_LEFT_FEETFIRST_CORONAL : slice_projection = MDC_CORONAL; break; default: slice_projection = MDC_TRANSAXIAL; } } return(slice_projection); } Int8 MdcGetNewPatSliceOrient(FILEINFO *cur, Int8 newproj) { Int8 pat_slice_orient=MDC_UNKNOWN; switch (cur->pat_slice_orient) { case MDC_SUPINE_HEADFIRST_TRANSAXIAL : case MDC_SUPINE_HEADFIRST_SAGITTAL : case MDC_SUPINE_HEADFIRST_CORONAL : switch (newproj) { case MDC_TRANSAXIAL : pat_slice_orient = MDC_SUPINE_HEADFIRST_TRANSAXIAL; break; case MDC_SAGITTAL : pat_slice_orient = MDC_SUPINE_HEADFIRST_SAGITTAL; break; case MDC_CORONAL : pat_slice_orient = MDC_SUPINE_HEADFIRST_CORONAL; break; } break; case MDC_PRONE_HEADFIRST_TRANSAXIAL : case MDC_PRONE_HEADFIRST_SAGITTAL : case MDC_PRONE_HEADFIRST_CORONAL : switch (newproj) { case MDC_TRANSAXIAL: pat_slice_orient = MDC_PRONE_HEADFIRST_TRANSAXIAL; break; case MDC_SAGITTAL : pat_slice_orient = MDC_PRONE_HEADFIRST_SAGITTAL; break; case MDC_CORONAL : pat_slice_orient = MDC_PRONE_HEADFIRST_CORONAL; break; } break; case MDC_SUPINE_FEETFIRST_TRANSAXIAL : case MDC_SUPINE_FEETFIRST_SAGITTAL : case MDC_SUPINE_FEETFIRST_CORONAL : switch (newproj) { case MDC_TRANSAXIAL : pat_slice_orient = MDC_SUPINE_FEETFIRST_TRANSAXIAL; break; case MDC_SAGITTAL : pat_slice_orient = MDC_SUPINE_FEETFIRST_SAGITTAL; break; case MDC_CORONAL : pat_slice_orient = MDC_SUPINE_FEETFIRST_CORONAL; break; } break; case MDC_PRONE_FEETFIRST_TRANSAXIAL : case MDC_PRONE_FEETFIRST_SAGITTAL : case MDC_PRONE_FEETFIRST_CORONAL : switch (newproj) { case MDC_TRANSAXIAL : pat_slice_orient = MDC_PRONE_FEETFIRST_TRANSAXIAL; break; case MDC_SAGITTAL : pat_slice_orient = MDC_PRONE_FEETFIRST_SAGITTAL; break; case MDC_CORONAL : pat_slice_orient = MDC_PRONE_FEETFIRST_CORONAL; break; } break; case MDC_DECUBITUS_RIGHT_HEADFIRST_TRANSAXIAL : case MDC_DECUBITUS_RIGHT_HEADFIRST_SAGITTAL : case MDC_DECUBITUS_RIGHT_HEADFIRST_CORONAL : switch (newproj) { case MDC_TRANSAXIAL : pat_slice_orient = MDC_DECUBITUS_RIGHT_HEADFIRST_TRANSAXIAL; break; case MDC_SAGITTAL : pat_slice_orient = MDC_DECUBITUS_RIGHT_HEADFIRST_SAGITTAL; break; case MDC_CORONAL : pat_slice_orient = MDC_DECUBITUS_RIGHT_HEADFIRST_CORONAL; break; } break; case MDC_DECUBITUS_LEFT_HEADFIRST_TRANSAXIAL : case MDC_DECUBITUS_LEFT_HEADFIRST_SAGITTAL : case MDC_DECUBITUS_LEFT_HEADFIRST_CORONAL : switch (newproj) { case MDC_TRANSAXIAL: pat_slice_orient = MDC_DECUBITUS_LEFT_HEADFIRST_TRANSAXIAL; break; case MDC_SAGITTAL : pat_slice_orient = MDC_DECUBITUS_LEFT_HEADFIRST_SAGITTAL; break; case MDC_CORONAL : pat_slice_orient = MDC_DECUBITUS_LEFT_HEADFIRST_CORONAL; break; } break; case MDC_DECUBITUS_RIGHT_FEETFIRST_TRANSAXIAL : case MDC_DECUBITUS_RIGHT_FEETFIRST_SAGITTAL : case MDC_DECUBITUS_RIGHT_FEETFIRST_CORONAL : switch (newproj) { case MDC_TRANSAXIAL : pat_slice_orient = MDC_DECUBITUS_RIGHT_FEETFIRST_TRANSAXIAL; break; case MDC_SAGITTAL : pat_slice_orient = MDC_DECUBITUS_RIGHT_FEETFIRST_SAGITTAL; break; case MDC_CORONAL : pat_slice_orient = MDC_DECUBITUS_RIGHT_FEETFIRST_CORONAL; break; } break; case MDC_DECUBITUS_LEFT_FEETFIRST_TRANSAXIAL : case MDC_DECUBITUS_LEFT_FEETFIRST_SAGITTAL : case MDC_DECUBITUS_LEFT_FEETFIRST_CORONAL : switch (newproj) { case MDC_TRANSAXIAL : pat_slice_orient = MDC_DECUBITUS_LEFT_FEETFIRST_TRANSAXIAL; break; case MDC_SAGITTAL : pat_slice_orient = MDC_DECUBITUS_LEFT_FEETFIRST_SAGITTAL; break; case MDC_CORONAL : pat_slice_orient = MDC_DECUBITUS_LEFT_FEETFIRST_CORONAL; break; } break; } return(pat_slice_orient); } char *MdcCheckReslice(FILEINFO *cur, Int8 newproj) { Int8 curproj; curproj = MdcGetSliceProjection(cur); /* some sanity checks before allowing reslicing */ if (cur->planar == MDC_YES) { strcpy(mdcbufr,"Planar study inappropriate"); return(mdcbufr); } /* don't fail in a batched job */ if (XMDC_GUI == MDC_YES) { if (newproj == curproj) { switch (curproj) { case MDC_TRANSAXIAL : sprintf(mdcbufr,"Already in XY - TRANSVERSE projection"); break; case MDC_SAGITTAL : sprintf(mdcbufr,"Already in YZ - SAGITTAL projection"); break; case MDC_CORONAL : sprintf(mdcbufr,"Already in XZ - CORONAL projection"); break; } return(mdcbufr); } } if (curproj == MDC_UNKNOWN) { strcpy(mdcbufr,"Current projection unknown"); return(mdcbufr); } if (cur->diff_type == MDC_YES) { strcpy(mdcbufr,"Identical pixel types required"); return(mdcbufr); } if (cur->diff_size == MDC_YES) { strcpy(mdcbufr,"Identical image sizes required"); return(mdcbufr); } if (cur->dim[3] <= 2) { strcpy(mdcbufr,"No volume detected"); return(mdcbufr); } if (cur->dim[3] <= MDC_VOLUME_REQUIRED) { strcpy(mdcbufr,"Volume too small"); return(mdcbufr); } if (cur->reconstructed == MDC_NO) { strcpy(mdcbufr,"Reconstructed data required"); return(mdcbufr); } return(NULL); } /* X, Y, Z: the true indices for x,y,z in our arrays */ /* OX, OY, OZ: old dim = new dim (do the reslice) */ /* DX, DY, DZ: new dim = old dim (get the sizes ) */ char *MdcResliceImages(FILEINFO *cur, Int8 newproj) { FILEINFO *new; IMG_DATA *newid, *curid; Uint32 nbytes, obytes, pixels, olength; Uint32 newX, newY, newZ, curX=0, curY=0, curZ=0, f, frames; Uint32 X=1, Y=2, Z=3, OX=X, OY=Y, OZ=Z, DX=X, DY=Y, DZ=Z; Uint8 *newp, *curp; Int8 curproj; double pixval; char *msg; curproj = MdcGetSliceProjection(cur); /* some sanity checks before doing reslice */ msg = MdcCheckReslice(cur,newproj); if (msg != NULL) return(msg); /* get temporary FILEINFO structure */ new = (FILEINFO *)malloc(sizeof(FILEINFO)); if (new == NULL) return("Couldn't malloc FILEINFO struct"); MdcCopyFI(new,cur,MDC_NO,MDC_YES); /* change orientation information */ new->pat_slice_orient = MdcGetNewPatSliceOrient(cur,newproj); strcpy(new->pat_orient,MdcGetStrPatOrient(new->pat_slice_orient)); /* prepare dimension mappings */ switch (newproj) { case MDC_TRANSAXIAL : switch (curproj) { case MDC_TRANSAXIAL : OX=X; OY=Y; OZ=Z; /* T -> T (===) */ DX=X; DY=Y; DZ=Z; break; case MDC_SAGITTAL : OX=Y; OY=Z; OZ=X; /* S -> T (sag) */ DX=Z; DY=X; DZ=Y; break; case MDC_CORONAL : OX=X; OY=Z; OZ=Y; /* C -> T (cor) */ DX=X; DY=Z; DZ=Y; break; } break; case MDC_SAGITTAL : switch (curproj) { case MDC_TRANSAXIAL : OX=Z; OY=X; OZ=Y; /* T -> S (sag) */ DX=Y; DY=Z; DZ=X; break; case MDC_SAGITTAL : OX=X; OY=Y; OZ=Z; /* S -> S (===) */ DX=X; DY=Y; DZ=Z; break; case MDC_CORONAL : OX=Z; OY=Y; OZ=X; /* C -> S (c2s) */ DX=Z; DY=Y; DZ=X; break; } break; case MDC_CORONAL : switch (curproj) { case MDC_TRANSAXIAL : OX=X; OY=Z; OZ=Y; /* T -> C (cor) */ DX=X; DY=Z; DZ=Y; break; case MDC_SAGITTAL : OX=Z; OY=Y; OZ=X; /* S -> C (s2c) */ DX=Z; DY=Y; DZ=X; break; case MDC_CORONAL : OX=X; OY=Y; OZ=Z; /* C -> C (===) */ DX=X; DY=Y; DZ=Z; break; } break; } /* first remove gaps between slices */ cur->pixdim[Z] = cur->image[0].slice_spacing; /* number of frames stay the same */ new->dim[0] = cur->dim[0]; for (frames=1, f=4; f<=cur->dim[0]; f++) { frames *= cur->dim[f]; new->dim[f] = cur->dim[f]; } new->pixdim[0] = cur->pixdim[0]; for (f=4; f<=cur->pixdim[0]; f++) new->pixdim[f] = cur->pixdim[f]; /* fill in new dimension values */ if (OX == 3) new->number = cur->dim[X] * frames; else if (OY == 3) new->number = cur->dim[Y] * frames; else if (OZ == 3) new->number = cur->dim[Z] * frames; new->dim[X] = cur->dim[DX]; new->pixdim[X] = cur->pixdim[DX]; new->dim[Y] = cur->dim[DY]; new->pixdim[Y] = cur->pixdim[DY]; new->dim[Z] = cur->dim[DZ]; new->pixdim[Z] = cur->pixdim[DZ]; new->mwidth = new->dim[X]; new->mheight = new->dim[Y]; /* handle pixel stuff */ if (MDC_QUANTIFY || MDC_CALIBRATE) { new->type = FLT32; new->bits = MdcType2Bits(new->type); }else{ new->type = cur->type; new->bits = cur->bits; } pixels = new->dim[X] * new->dim[Y]; obytes = MdcType2Bytes(cur->type); nbytes = MdcType2Bytes(new->type); /* get new IMG_DATA structs */ if (!MdcGetStructID(new,new->number)) { MdcCleanUpFI(new); return("Couldn't malloc IMG_DATA structs"); } /* reslice images and fill in structures */ olength = cur->dim[X]; for (f=0; fdim[Z]; newZ++) { newid = &new->image[newZ + (f * new->dim[Z])]; newid->buf = MdcGetImgBuffer(pixels * nbytes); if (newid->buf == NULL) { MdcCleanUpFI(new); return("Couldn't malloc image buffer"); } newp = newid->buf; newid->width = new->mwidth; newid->height= new->mheight; newid->bits = new->bits; newid->type = new->type; newid->pixel_xsize = new->pixdim[X]; newid->pixel_ysize = new->pixdim[Y]; newid->slice_width = new->pixdim[Z]; newid->slice_spacing= new->pixdim[Z]; MdcFillImgPos(new,newZ,newZ,0.); MdcFillImgOrient(new,newZ); for (newY=0; newYdim[Y]; newY++) { for (newX=0; newXdim[X]; newX++) { /* X-mapping */ if (OX == X) curX = newX; else if (OX == Y) curX = newY; else if (OX == Z) curX = newZ; /* Y-mapping */ if (OY == X) curY = newX; else if (OY == Y) curY = newY; else if (OY == Z) curY = newZ; /* Z-mapping */ if (OZ == X) curZ = newX; else if (OZ == Y) curZ = newY; else if (OZ == Z) curZ = newZ; curid = &cur->image[curZ + (f * cur->dim[Z])]; curp = curid->buf + (((curY * olength) + curX) * obytes); if (MDC_QUANTIFY || MDC_CALIBRATE) { pixval = MdcGetDoublePixel(curp,curid->type); pixval *= (double)curid->rescale_slope; pixval += (double)curid->rescale_intercept; MdcPutDoublePixel(newp,pixval,newid->type); }else{ memcpy(newp,curp,nbytes); } newp += nbytes; } } } if (cur->acquisition_type == MDC_ACQUISITION_GATED || cur->acquisition_type == MDC_ACQUISITION_GSPECT) { /* alter a gated parameter to keep things in line */ /* through different reslices: based on HeartRate */ /* since this depends on number of images per frame */ if (cur->gdata != NULL && new->gdata != NULL) { new->gdata[0].time_per_proj = cur->gdata[0].time_per_proj; new->gdata[0].time_per_proj *= (float)cur->dim[Z]; new->gdata[0].time_per_proj /= (float)new->dim[Z]; } } if (cur->acquisition_type == MDC_ACQUISITION_TOMO || cur->acquisition_type == MDC_ACQUISITION_DYNAMIC) { /* Just fix nr_of_slices for each DYNAMIC_DATA struct. */ /* All other entries can be preserved for tomo study. */ for (f=0; fdynnr; f++) new->dyndata[f].nr_of_slices = new->dim[3]; } /* set new slice projection */ new->slice_projection = newproj; /* check integrity */ if ((msg = MdcImagesPixelFiddle(new)) != NULL) { MdcCleanUpFI(new); MdcFree(new) return(msg); } /* remove cur */ MdcCleanUpFI(cur); /* copy new -> cur */ MdcCopyFI(cur,new,MDC_NO,MDC_YES); /* just rehang image pointer */ cur->number= new->number; cur->image = new->image; /* and mask new image pointer */ new->number = 0; new->image = NULL; /* now safely remove new */ MdcCleanUpFI(new); MdcFree(new); return(NULL); }