heart.c File Reference

Functions for simple simulation of image of heart. More...

#include "libtpcidi.h"

Go to the source code of this file.

Functions
int	imgSimulateRing (IMG *img, int fi, int zi, double cx, double cy, double r1, double r2, double vr, double vi, double vo, int verbose)
int	simMyocDiameterCurve (DFT *dft, double t1, double t2, double hbr, double maxdiam, double mindiam)
int	imgSimulateSphere (IMG *img, int fi, double cx, double cy, double cz, double r1, double r2, double vr, double vi, double vo, int verbose)

Detailed Description

Functions for simple simulation of image of heart.

Author

Vesa Oikonen

Definition in file heart.c.

Function Documentation

imgSimulateRing()

int imgSimulateRing	(	IMG *	img,
		int	fi,
		int	zi,
		double	cx,
		double	cy,
		double	r1,
		double	r2,
		double	vr,
		double	vi,
		double	vo,
		int	verbose )

Simulate an image of ring with specified inner and outer radius and activity in the ring and inside and outside of it. The applied method is only approximate at pixel borders (pixel is divided into 5x5 subpixels). PET resolution effects are not simulated.

See also

imgRingMask, imgSimulateSphere, imgGaussianFIRFilter, idiSimulateTubeVol

Returns

Returns 0 if successful.

Parameters

img	Pointer to allocated static or dynamic image; image must contain pixel sizes and dimensions; values are added to any existing pixel values, thus you may need to set pixel values to zero before calling this function.
fi	Frame index [0..dimt-1].
zi	Plane index [0..dimz-1].
cx	X distance of circle centre (mm) from the upper left corner of the image.
cy	Y distance of circle centre (mm) from the upper left corner of the image.
r1	Inner radius of circle (mm).
r2	Outer radius of circle (mm).
vr	Ring value; this value is added to each pixel value that fits inside the radius's; the pixels that are partially inside the ring will get fraction of the value.
vi	Inside value; this value is added to each pixel value that fits inside the inner radius; the pixels that are partially inside the radius will get fraction of the value.
vo	Outside value; this value is added to each pixel value that fits outside the outer radius; the pixels that are partially outside the radius will get fraction of the value.
verbose	Verbose level; set to <=0 to prevent all prints to stdout.

Definition at line 17 of file heart.c.

  {
  int nr, ni, no, xi, yi, i, j;
  double dx[5], dy[5], v, d;
  double ri2, ro2;
 
  if(verbose>0)
     printf("%s(img, %d, %d, %g, %g, %g, %g, %g, %g, %g, %d)\n",
            __func__, fi, zi, cx, cy, r1, r2, vr, vi, vo, verbose);
  if(r1>=r2 || r1<0.0) {
    if(verbose>0) fprintf(stderr, "Error: invalid radius.\n");
    return(1);
  }
 
  if(img->status<IMG_STATUS_OCCUPIED) return(1);
  if(zi<0 || zi>=img->dimz) return(2);
  if(img->sizey<=0.0 || img->sizex<=0.0) return(3);
  if(img->sizey!=img->sizex) return(4);
  if(img->dimt<1 || fi>=img->dimt) return(5);
  if(r1<0.0 || r2<r1) return(6);
 
  /* Compare radius ^2 to sum of squared distances instead of square roots */
  ri2=r1*r1; ro2=r2*r2;
  for(yi=0; yi<img->dimy; yi++) {
    dy[0]=(0.1+(double)yi)*img->sizey - cy;
    for(j=1; j<5; j++) dy[j]=dy[j-1]+0.2*img->sizey;
    for(xi=0; xi<img->dimx; xi++) {
      dx[0]=(0.1+(double)xi)*img->sizex - cx;
      for(i=1; i<5; i++) dx[i]=dx[i-1]+0.2*img->sizex;
      nr=ni=no=0;
      for(i=0; i<5; i++) for(j=0; j<5; j++) {
        d=dx[i]*dx[i]+dy[j]*dy[j];
        if(d<ri2) ni++; else if(d<ro2) nr++; else no++;
      }
      v=(double)nr*vr/25.0;
      v+=(double)ni*vi/25.0;
      v+=(double)no*vo/25.0;
      img->m[zi][yi][xi][fi]+=v;
    }
  }
  return(0);
}

imgSimulateSphere()

int imgSimulateSphere	(	IMG *	img,
		int	fi,
		double	cx,
		double	cy,
		double	cz,
		double	r1,
		double	r2,
		double	vr,
		double	vi,
		double	vo,
		int	verbose )

Simulate a 3D image of circle with specified inner and outer radius and activity in the circle and inside and outside of it.

The applied method is only approximate at pixel borders (pixel is divided into 5x5 subpixels). PET resolution effects are not simulated.

See also

imgRingMask, imgCircleMask, imgGaussianFIRFilter, idiSimulateTubeVol

Returns

Returns 0 if successful.

Parameters

img	Pointer to allocated static or dynamic image; image must contain pixel sizes and dimensions; values are added to any existing pixel values, thus you may need to set pixel values to zero before calling this function.
fi	Frame index [0..dimt-1].
cx	X distance of circle centre (mm) from the upper left corner of the image.
cy	Y distance of circle centre (mm) from the upper left corner of the image.
cz	Z distance of circle centre (mm) from the upper left corner of the image.
r1	Inner radius of circle (mm).
r2	Outer radius of circle (mm).
vr	Sphere wall value; this value is added to each pixel value that fits inside the radius's; the pixels that are partially inside the wall will get fraction of the value.
vi	Inside value; this value is added to each pixel value that fits inside the inner radius; the pixels that are partially inside the radius will get fraction of the value.
vo	Outside value; this value is added to each pixel value that fits outside the outer radius; the pixels that are partially outside the radius will get fraction of the value.
verbose	Verbose level; set to <=0 to prevent all prints to stdout.

Definition at line 154 of file heart.c.

  {
  int nr, ni, no, xi, yi, zi, i, j, k;
  double dx[5], dy[5], dz[5], v, d;
  double ri2, ro2;
 
  if(verbose>0)
     printf("%s(img, %d, %g, %g, %g, %g, %g, %g, %g, %g, %d)\n",
            __func__, fi, cx, cy, cz, r1, r2, vr, vi, vo, verbose);
  if(r1>=r2 || r1<0.0) {
    if(verbose>0) fprintf(stderr, "Error: invalid radius.\n");
    return(1);
  }
 
  if(img->status<IMG_STATUS_OCCUPIED) return(1);
  if(img->sizez<=0.0 || img->sizey<=0.0 || img->sizex<=0.0) return(2);
  if(img->sizey!=img->sizex || img->sizez!=img->sizex) return(3);
  if(img->dimt<1 || fi>=img->dimt) return(4);
  if(r1<0.0 || r2<r1) return(5);
 
  /* Compare radius ^2 to sum of squared distances instead of square roots */
  ri2=r1*r1; ro2=r2*r2;
  for(zi=0; zi<img->dimz; zi++) {
    dz[0]=(0.1+(double)zi)*img->sizez - cz;
    for(k=1; k<5; k++) dz[k]=dz[k-1]+0.2*img->sizez;
    for(yi=0; yi<img->dimy; yi++) {
      dy[0]=(0.1+(double)yi)*img->sizey - cy;
      for(j=1; j<5; j++) dy[j]=dy[j-1]+0.2*img->sizey;
      for(xi=0; xi<img->dimx; xi++) {
        dx[0]=(0.1+(double)xi)*img->sizex - cx;
        for(i=1; i<5; i++) dx[i]=dx[i-1]+0.2*img->sizex;
        nr=ni=no=0;
        for(i=0; i<5; i++) for(j=0; j<5; j++) for(k=0; k<5; k++) {
          d=dx[i]*dx[i]+dy[j]*dy[j]+dz[k]*dz[k];
          if(d<ri2) ni++; else if(d<ro2) nr++; else no++;
        }
        v=(double)nr*vr/125.0;
        v+=(double)ni*vi/125.0;
        v+=(double)no*vo/125.0;
        img->m[zi][yi][xi][fi]+=v;
      }
    }
  }
  return(0);
}

simMyocDiameterCurve()

int simMyocDiameterCurve	(	DFT *	dft,
		double	t1,
		double	t2,
		double	hbr,
		double	maxdiam,
		double	mindiam )

Calculate the inner LV cavity diameter as a function of time for simulations.

10 samples per heart beat will be calculated. Fractional increase/decrease rate is fixed and coded in here.

See also

imgSimulateRing, imgRingMask, idiSimulateTubeVol

Returns

Returns 0 if successful.

Parameters

dft	Pointer to initiated DFT struct where diameter curve will be written.
t1	Start time (s).
t2	Stop time (s).
hbr	Heart rate (beats/min).
maxdiam	Maximum inner diameter (mm).
mindiam	Minimum inner diameter (mm).

Definition at line 96 of file heart.c.

  {
  int n, i, ret;
  double s;
  double rfract[10]={1.0, 0.6, 0.25, 0.0, 0.1, 0.25, 0.6, 0.8, 0.9, 0.95};
 
  /* check input */
  if(dft==NULL) return(1);
  if(t2<t1) return(1);
  if(maxdiam<mindiam) return(1);
 
  /* delete any previous curve contents */
  dftEmpty(dft);
 
  /* convert heart rate to per sec */
  if(hbr>0.0) hbr/=60.0; else hbr=0.0;
 
  /* calculate the nr of samples */
  if(hbr>0.0 && (t2-t1)>0.0) n=10.0*(t2-t1)/hbr; else n=1;
 
  /* Set up curve struct */
  ret=dftSetmem(dft, n, 1); if(ret!=0) return(3);
  dft->frameNr=n; dft->voiNr=1; dft->timetype=DFT_TIME_STARTEND;
  dft->timeunit=TUNIT_SEC;
 
  /* Fill curve */
  dft->x1[0]=t1; if(hbr>0.0) s=0.1*hbr; else s=t2-t1; 
  for(i=n=0; i<dft->frameNr; i++) {
    if(i>0) dft->x1[i]=dft->x2[i-1];
    dft->x2[i]=dft->x1[i]+s; dft->x[i]=0.5*(dft->x1[i]+dft->x2[i]);
    dft->voi[0].y[i]=mindiam+rfract[n]*(maxdiam-mindiam); n++; if(n==10) n=0;
  }
 
  return 0;
}