peak.c File Reference

Functions for finding peak voxel in PET image. More...

#include "libtpcidi.h"

Go to the source code of this file.

Functions
int	imgGetConcWeightedPeakPos (IMG *img, float thrs, IMG_PIXEL *pos, int verbose)

Detailed Description

Functions for finding peak voxel in PET image.

Author

Vesa Oikonen

Definition in file peak.c.

Function Documentation

imgGetConcWeightedPeakPos()

int imgGetConcWeightedPeakPos	(	IMG *	img,
		float	thrs,
		IMG_PIXEL *	pos,
		int	verbose )

Calculates concentration weighted peak position from a sum image.

See also

imgMaskPixelTACs, imgCircleMask, idiSimulateTubeVol

Returns

Returns 0 if successful.

Parameters

img	Pointer to allocated image from which the peak position is searched. If IMG struct contains more than one times frame, then only the first one is used.
thrs	Threshold fraction (0<thrs<1) between image min and max.
pos	Pointer to pixel position structure.
verbose	Verbose level; set to <=0 to prevent all prints to stdout and stderr.

Definition at line 14 of file peak.c.

  {
  if(verbose>0) printf("%s(img, %g, ...)\n", __func__, thrs);
 
  int xi, yi, zi, fi=0;
  float d, fmax, fmin, thrslev;
  float px, py, pz, fsum;
 
  if(img->status<IMG_STATUS_OCCUPIED) return(1);
  if(img->dimz<1 || img->dimx<2 || img->dimy<2) return(2);
  if(img->dimt<1) return(3);
  if(thrs<=0.0 || thrs>=1.0) return(4);
  if(pos==NULL) return(5);
 
 
  /* Search min and max value inside the image volume */
  zi=yi=xi=0; fmax=fmin=img->m[zi][yi][xi][fi];
  for(zi=0; zi<img->dimz; zi++)
    for(yi=0; yi<img->dimy; yi++)
      for(xi=0; xi<img->dimx; xi++) {
        if(img->m[zi][yi][xi][fi] > fmax) fmax=img->m[zi][yi][xi][fi];
        else if(img->m[zi][yi][xi][fi] < fmin) fmin=img->m[zi][yi][xi][fi];
      }
  if(verbose>1) printf("fmin := %g\nfmax := %g\n", fmin, fmax);
 
  /* Calculate threshold level from min and max */
  d=fmax-fmin; if(d<=0.0) return(10);
  thrslev=d*thrs+fmin; if(verbose>1) printf("thrslev := %g\n", thrslev);
 
  /* Calculate concentration weighted peak position */
  px=py=pz=0.0; fsum=0.0;
  for(zi=0; zi<img->dimz; zi++)
    for(yi=0; yi<img->dimy; yi++)
      for(xi=0; xi<img->dimx; xi++) {
        d=img->m[zi][yi][xi][fi]-thrslev;
        if(d>0.0) {
          px+=(float)(xi+1)*d;
          py+=(float)(yi+1)*d;
          pz+=(float)(zi+1)*d;
          fsum+=d;
        }
      }
  if(fsum<=0.0) return 5;
  px/=fsum; py/=fsum; pz/=fsum;
  if(verbose>1) printf("peak_pos := (%g,%g,%g)\n", px, py, pz);
  pos->x=roundf(px); pos->y=roundf(py); pos->z=roundf(pz); pos->f=1;
 
  return(0);
}